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

Last change on this file since 534 was 494, checked in by raasch, 15 years ago

last commit documented; configuration example file for netcdf4 added

  • Property svn:keywords set to Id
File size: 76.3 KB
Line 
1 SUBROUTINE header
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: header.f90 494 2010-03-01 12:23:32Z heinze $
11!
12! 493 2010-03-01 08:30:24Z raasch
13! NetCDF data output format extendend for NetCDF4/HDF5
14!
15! 449 2010-02-02 11:23:59Z raasch
16! +large scale vertical motion (subsidence/ascent)
17! Bugfix: index problem concerning gradient_level indices removed
18!
19! 410 2009-12-04 17:05:40Z letzel
20! Masked data output: + dt_domask, mask_01~20_x|y|z, mask_01~20_x|y|z_loop,
21! mask_scale|_x|y|z, masks, skip_time_domask
22!
23! 346 2009-07-06 10:13:41Z raasch
24! initializing_actions='read_data_for_recycling' renamed to 'cyclic_fill'
25! Coupling with independent precursor runs.
26! Output of messages replaced by message handling routine.
27! Output of several additional dvr parameters
28! +canyon_height, canyon_width_x, canyon_width_y, canyon_wall_left,
29! canyon_wall_south, conserve_volume_flow_mode, dp_external, dp_level_b,
30! dp_smooth, dpdxy, u_bulk, v_bulk
31! topography_grid_convention moved from user_header
32! small bugfix concerning 3d 64bit netcdf output format
33!
34! 206 2008-10-13 14:59:11Z raasch
35! Bugfix: error in zu index in case of section_xy = -1
36!
37! 198 2008-09-17 08:55:28Z raasch
38! Format adjustments allowing output of larger revision numbers
39!
40! 197 2008-09-16 15:29:03Z raasch
41! allow 100 spectra levels instead of 10 for consistency with
42! define_netcdf_header,
43! bugfix in the output of the characteristic levels of potential temperature,
44! geostrophic wind, scalar concentration, humidity and leaf area density,
45! output of turbulence recycling informations
46!
47! 138 2007-11-28 10:03:58Z letzel
48! Allow new case bc_uv_t = 'dirichlet_0' for channel flow.
49! Allow two instead of one digit to specify isosurface and slicer variables.
50! Output of sorting frequency of particles
51!
52! 108 2007-08-24 15:10:38Z letzel
53! Output of informations for coupled model runs (boundary conditions etc.)
54! + output of momentumfluxes at the top boundary
55! Rayleigh damping for ocean, e_init
56!
57! 97 2007-06-21 08:23:15Z raasch
58! Adjustments for the ocean version.
59! use_pt_reference renamed use_reference
60!
61! 87 2007-05-22 15:46:47Z raasch
62! Bugfix: output of use_upstream_for_tke
63!
64! 82 2007-04-16 15:40:52Z raasch
65! Preprocessor strings for different linux clusters changed to "lc",
66! routine local_flush is used for buffer flushing
67!
68! 76 2007-03-29 00:58:32Z raasch
69! Output of netcdf_64bit_3d, particles-package is now part of the default code,
70! output of the loop optimization method, moisture renamed humidity,
71! output of subversion revision number
72!
73! 19 2007-02-23 04:53:48Z raasch
74! Output of scalar flux applied at top boundary
75!
76! RCS Log replace by Id keyword, revision history cleaned up
77!
78! Revision 1.63  2006/08/22 13:53:13  raasch
79! Output of dz_max
80!
81! Revision 1.1  1997/08/11 06:17:20  raasch
82! Initial revision
83!
84!
85! Description:
86! ------------
87! Writing a header with all important informations about the actual run.
88! This subroutine is called three times, two times at the beginning
89! (writing information on files RUN_CONTROL and HEADER) and one time at the
90! end of the run, then writing additional information about CPU-usage on file
91! header.
92!-----------------------------------------------------------------------------!
93
94    USE arrays_3d
95    USE control_parameters
96    USE cloud_parameters
97    USE cpulog
98    USE dvrp_variables
99    USE grid_variables
100    USE indices
101    USE model_1d
102    USE particle_attributes
103    USE pegrid
104    USE subsidence_mod
105    USE spectrum
106
107    IMPLICIT NONE
108
109    CHARACTER (LEN=1)  ::  prec
110    CHARACTER (LEN=2)  ::  do2d_mode
111    CHARACTER (LEN=5)  ::  section_chr
112    CHARACTER (LEN=9)  ::  time_to_string
113    CHARACTER (LEN=10) ::  coor_chr, host_chr
114    CHARACTER (LEN=16) ::  begin_chr
115    CHARACTER (LEN=23) ::  ver_rev
116    CHARACTER (LEN=40) ::  output_format
117    CHARACTER (LEN=70) ::  char1, char2, dopr_chr, &
118                           do2d_xy, do2d_xz, do2d_yz, do3d_chr, &
119                           domask_chr, run_classification
120    CHARACTER (LEN=86) ::  coordinates, gradients, learde, slices,  &
121                           temperatures, ugcomponent, vgcomponent
122    CHARACTER (LEN=85) ::  roben, runten
123
124    CHARACTER (LEN=1), DIMENSION(1:3) ::  dir = (/ 'x', 'y', 'z' /)
125
126    INTEGER ::  av, bh, blx, bly, bxl, bxr, byn, bys, ch, count, cwx, cwy,  &
127         cxl, cxr, cyn, cys, dim, i, ihost, io, j, l, ll, m, mpi_type
128    REAL    ::  cpuseconds_per_simulated_second
129
130!
131!-- Open the output file. At the end of the simulation, output is directed
132!-- to unit 19.
133    IF ( ( runnr == 0 .OR. force_print_header )  .AND. &
134         .NOT. simulated_time_at_begin /= simulated_time )  THEN
135       io = 15   !  header output on file RUN_CONTROL
136    ELSE
137       io = 19   !  header output on file HEADER
138    ENDIF
139    CALL check_open( io )
140
141!
142!-- At the end of the run, output file (HEADER) will be rewritten with
143!-- new informations
144    IF ( io == 19 .AND. simulated_time_at_begin /= simulated_time ) REWIND( 19 )
145
146!
147!-- Determine kind of model run
148    IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
149       run_classification = '3D - restart run'
150    ELSEIF ( TRIM( initializing_actions ) == 'cyclic_fill' )  THEN
151       run_classification = '3D - run with cyclic fill of 3D - prerun data'
152    ELSEIF ( INDEX( initializing_actions, 'set_constant_profiles' ) /= 0 )  THEN
153       run_classification = '3D - run without 1D - prerun'
154    ELSEIF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
155       run_classification = '3D - run with 1D - prerun'
156    ELSEIF ( INDEX( initializing_actions, 'by_user' ) /=0 )  THEN
157       run_classification = '3D - run initialized by user'
158    ELSE
159       message_string = ' unknown action(s): ' // TRIM( initializing_actions )
160       CALL message( 'header', 'PA0191', 0, 0, 0, 6, 0 )
161    ENDIF
162    IF ( ocean )  THEN
163       run_classification = 'ocean - ' // run_classification
164    ELSE
165       run_classification = 'atmosphere - ' // run_classification
166    ENDIF
167
168!
169!-- Run-identification, date, time, host
170    host_chr = host(1:10)
171    ver_rev = TRIM( version ) // '  ' // TRIM( revision )
172    WRITE ( io, 100 )  ver_rev, TRIM( run_classification )
173    IF ( TRIM( coupling_mode ) /= 'uncoupled' )  THEN
174#if defined( __mpi2 )
175       mpi_type = 2
176#else
177       mpi_type = 1
178#endif
179       WRITE ( io, 101 )  mpi_type, coupling_mode
180    ENDIF
181    WRITE ( io, 102 )  run_date, run_identifier, run_time, runnr, &
182                       ADJUSTR( host_chr )
183#if defined( __parallel )
184    IF ( npex == -1  .AND.  pdims(2) /= 1 )  THEN
185       char1 = 'calculated'
186    ELSEIF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.  &
187               host(1:2) == 'lc' )  .AND.                          &
188             npex == -1  .AND.  pdims(2) == 1 )  THEN
189       char1 = 'forced'
190    ELSE
191       char1 = 'predefined'
192    ENDIF
193    IF ( threads_per_task == 1 )  THEN
194       WRITE ( io, 103 )  numprocs, pdims(1), pdims(2), TRIM( char1 )
195    ELSE
196       WRITE ( io, 104 )  numprocs*threads_per_task, numprocs, &
197                          threads_per_task, pdims(1), pdims(2), TRIM( char1 )
198    ENDIF
199    IF ( ( host(1:3) == 'ibm'  .OR.  host(1:3) == 'nec'  .OR.    &
200           host(1:2) == 'lc'   .OR.  host(1:3) == 'dec' )  .AND. &
201         npex == -1  .AND.  pdims(2) == 1 )                      &
202    THEN
203       WRITE ( io, 106 )
204    ELSEIF ( pdims(2) == 1 )  THEN
205       WRITE ( io, 107 )  'x'
206    ELSEIF ( pdims(1) == 1 )  THEN
207       WRITE ( io, 107 )  'y'
208    ENDIF
209    IF ( use_seperate_pe_for_dvrp_output )  WRITE ( io, 105 )
210#endif
211    WRITE ( io, 99 )
212
213!
214!-- Numerical schemes
215    WRITE ( io, 110 )
216    IF ( psolver(1:7) == 'poisfft' )  THEN
217       WRITE ( io, 111 )  TRIM( fft_method )
218       IF ( psolver == 'poisfft_hybrid' )  WRITE ( io, 138 )
219    ELSEIF ( psolver == 'sor' )  THEN
220       WRITE ( io, 112 )  nsor_ini, nsor, omega_sor
221    ELSEIF ( psolver == 'multigrid' )  THEN
222       WRITE ( io, 135 )  cycle_mg, maximum_grid_level, ngsrb
223       IF ( mg_cycles == -1 )  THEN
224          WRITE ( io, 140 )  residual_limit
225       ELSE
226          WRITE ( io, 141 )  mg_cycles
227       ENDIF
228       IF ( mg_switch_to_pe0_level == 0 )  THEN
229          WRITE ( io, 136 )  nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
230                             nzt_mg(1)
231       ELSEIF (  mg_switch_to_pe0_level /= -1 )  THEN
232          WRITE ( io, 137 )  mg_switch_to_pe0_level,            &
233                             mg_loc_ind(2,0)-mg_loc_ind(1,0)+1, &
234                             mg_loc_ind(4,0)-mg_loc_ind(3,0)+1, &
235                             nzt_mg(mg_switch_to_pe0_level),    &
236                             nxr_mg(1)-nxl_mg(1)+1, nyn_mg(1)-nys_mg(1)+1, &
237                             nzt_mg(1)
238       ENDIF
239    ENDIF
240    IF ( call_psolver_at_all_substeps  .AND. timestep_scheme(1:5) == 'runge' ) &
241    THEN
242       WRITE ( io, 142 )
243    ENDIF
244
245    IF ( momentum_advec == 'pw-scheme' )  THEN
246       WRITE ( io, 113 )
247    ELSE
248       WRITE ( io, 114 )
249       IF ( cut_spline_overshoot )  WRITE ( io, 124 )
250       IF ( overshoot_limit_u /= 0.0  .OR.  overshoot_limit_v /= 0.0  .OR. &
251            overshoot_limit_w /= 0.0 )  THEN
252          WRITE ( io, 127 )  overshoot_limit_u, overshoot_limit_v, &
253                             overshoot_limit_w
254       ENDIF
255       IF ( ups_limit_u /= 0.0  .OR.  ups_limit_v /= 0.0  .OR. &
256            ups_limit_w /= 0.0 )                               &
257       THEN
258          WRITE ( io, 125 )  ups_limit_u, ups_limit_v, ups_limit_w
259       ENDIF
260       IF ( long_filter_factor /= 0.0 )  WRITE ( io, 115 )  long_filter_factor
261    ENDIF
262    IF ( scalar_advec == 'pw-scheme' )  THEN
263       WRITE ( io, 116 )
264    ELSEIF ( scalar_advec == 'ups-scheme' )  THEN
265       WRITE ( io, 117 )
266       IF ( cut_spline_overshoot )  WRITE ( io, 124 )
267       IF ( overshoot_limit_e /= 0.0  .OR.  overshoot_limit_pt /= 0.0 )  THEN
268          WRITE ( io, 128 )  overshoot_limit_e, overshoot_limit_pt
269       ENDIF
270       IF ( ups_limit_e /= 0.0  .OR.  ups_limit_pt /= 0.0 )  THEN
271          WRITE ( io, 126 )  ups_limit_e, ups_limit_pt
272       ENDIF
273    ELSE
274       WRITE ( io, 118 )
275    ENDIF
276
277    WRITE ( io, 139 )  TRIM( loop_optimization )
278
279    IF ( galilei_transformation )  THEN
280       IF ( use_ug_for_galilei_tr )  THEN
281          char1 = 'geostrophic wind'
282       ELSE
283          char1 = 'mean wind in model domain'
284       ENDIF
285       IF ( simulated_time_at_begin == simulated_time )  THEN
286          char2 = 'at the start of the run'
287       ELSE
288          char2 = 'at the end of the run'
289       ENDIF
290       WRITE ( io, 119 )  TRIM( char1 ), TRIM( char2 ), &
291                          advected_distance_x/1000.0, advected_distance_y/1000.0
292    ENDIF
293    IF ( timestep_scheme == 'leapfrog' )  THEN
294       WRITE ( io, 120 )
295    ELSEIF ( timestep_scheme == 'leapfrog+euler' )  THEN
296       WRITE ( io, 121 )
297    ELSE
298       WRITE ( io, 122 )  timestep_scheme
299    ENDIF
300    IF ( use_upstream_for_tke )  WRITE ( io, 143 )
301    IF ( rayleigh_damping_factor /= 0.0 )  THEN
302       IF ( .NOT. ocean )  THEN
303          WRITE ( io, 123 )  'above', rayleigh_damping_height, &
304               rayleigh_damping_factor
305       ELSE
306          WRITE ( io, 123 )  'below', rayleigh_damping_height, &
307               rayleigh_damping_factor
308       ENDIF
309    ENDIF
310    IF ( humidity )  THEN
311       IF ( .NOT. cloud_physics )  THEN
312          WRITE ( io, 129 )
313       ELSE
314          WRITE ( io, 130 )
315          WRITE ( io, 131 )
316          IF ( radiation )      WRITE ( io, 132 )
317          IF ( precipitation )  WRITE ( io, 133 )
318       ENDIF
319    ENDIF
320    IF ( passive_scalar )  WRITE ( io, 134 )
321    IF ( conserve_volume_flow )  THEN
322       WRITE ( io, 150 )  conserve_volume_flow_mode
323       IF ( TRIM( conserve_volume_flow_mode ) == 'bulk_velocity' )  THEN
324          WRITE ( io, 151 )  u_bulk, v_bulk
325       ENDIF
326    ELSEIF ( dp_external )  THEN
327       IF ( dp_smooth )  THEN
328          WRITE ( io, 152 )  dpdxy, dp_level_b, ', vertically smoothed.'
329       ELSE
330          WRITE ( io, 152 )  dpdxy, dp_level_b, '.'
331       ENDIF
332    ENDIF
333    IF ( large_scale_subsidence )  THEN
334        WRITE ( io, 153 )
335        WRITE ( io, 154 )
336    ENDIF
337    WRITE ( io, 99 )
338
339!
340!-- Runtime and timestep informations
341    WRITE ( io, 200 )
342    IF ( .NOT. dt_fixed )  THEN
343       WRITE ( io, 201 )  dt_max, cfl_factor
344    ELSE
345       WRITE ( io, 202 )  dt
346    ENDIF
347    WRITE ( io, 203 )  simulated_time_at_begin, end_time
348
349    IF ( time_restart /= 9999999.9  .AND. &
350         simulated_time_at_begin == simulated_time )  THEN
351       IF ( dt_restart == 9999999.9 )  THEN
352          WRITE ( io, 204 )  ' Restart at:       ',time_restart
353       ELSE
354          WRITE ( io, 205 )  ' Restart at:       ',time_restart, dt_restart
355       ENDIF
356    ENDIF
357
358    IF ( simulated_time_at_begin /= simulated_time )  THEN
359       i = MAX ( log_point_s(10)%counts, 1 )
360       IF ( ( simulated_time - simulated_time_at_begin ) == 0.0 )  THEN
361          cpuseconds_per_simulated_second = 0.0
362       ELSE
363          cpuseconds_per_simulated_second = log_point_s(10)%sum / &
364                                            ( simulated_time -    &
365                                              simulated_time_at_begin )
366       ENDIF
367       WRITE ( io, 206 )  simulated_time, log_point_s(10)%sum, &
368                          log_point_s(10)%sum / REAL( i ),     &
369                          cpuseconds_per_simulated_second
370       IF ( time_restart /= 9999999.9  .AND.  time_restart < end_time )  THEN
371          IF ( dt_restart == 9999999.9 )  THEN
372             WRITE ( io, 204 )  ' Next restart at:  ',time_restart
373          ELSE
374             WRITE ( io, 205 )  ' Next restart at:  ',time_restart, dt_restart
375          ENDIF
376       ENDIF
377    ENDIF
378
379!
380!-- Start time for coupled runs, if independent precursor runs for atmosphere
381!-- and ocean are used. In this case, coupling_start_time defines the time
382!-- when the coupling is switched on.
383    IF ( coupling_start_time /= 0.0 )  THEN
384       IF ( coupling_start_time >= simulated_time_at_begin )  THEN
385          char1 = 'Precursor run for a coupled atmosphere-ocean run'
386       ELSE
387          char1 = 'Coupled atmosphere-ocean run following independent ' // &
388                  'precursor runs'
389       ENDIF
390       WRITE ( io, 207 )  char1, coupling_start_time
391    ENDIF
392
393!
394!-- Computational grid
395    IF ( .NOT. ocean )  THEN
396       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(nzt+1)
397       IF ( dz_stretch_level_index < nzt+1 )  THEN
398          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
399                             dz_stretch_factor, dz_max
400       ENDIF
401    ELSE
402       WRITE ( io, 250 )  dx, dy, dz, (nx+1)*dx, (ny+1)*dy, zu(0)
403       IF ( dz_stretch_level_index > 0 )  THEN
404          WRITE ( io, 252 )  dz_stretch_level, dz_stretch_level_index, &
405                             dz_stretch_factor, dz_max
406       ENDIF
407    ENDIF
408    WRITE ( io, 254 )  nx, ny, nzt+1, MIN( nnx, nx+1 ), MIN( nny, ny+1 ), &
409                       MIN( nnz+2, nzt+2 )
410    IF ( numprocs > 1 )  THEN
411       IF ( nxa == nx  .AND.  nya == ny  .AND.  nza == nz )  THEN
412          WRITE ( io, 255 )
413       ELSE
414          WRITE ( io, 256 )  nnx-(nxa-nx), nny-(nya-ny), nzt+2
415       ENDIF
416    ENDIF
417    IF ( sloping_surface )  WRITE ( io, 260 )  alpha_surface
418
419!
420!-- Topography
421    WRITE ( io, 270 )  topography
422    SELECT CASE ( TRIM( topography ) )
423
424       CASE ( 'flat' )
425          ! no actions necessary
426
427       CASE ( 'single_building' )
428          blx = INT( building_length_x / dx )
429          bly = INT( building_length_y / dy )
430          bh  = INT( building_height / dz )
431
432          IF ( building_wall_left == 9999999.9 )  THEN
433             building_wall_left = ( nx + 1 - blx ) / 2 * dx
434          ENDIF
435          bxl = INT ( building_wall_left / dx + 0.5 )
436          bxr = bxl + blx
437
438          IF ( building_wall_south == 9999999.9 )  THEN
439             building_wall_south = ( ny + 1 - bly ) / 2 * dy
440          ENDIF
441          bys = INT ( building_wall_south / dy + 0.5 )
442          byn = bys + bly
443
444          WRITE ( io, 271 )  building_length_x, building_length_y, &
445                             building_height, bxl, bxr, bys, byn
446
447       CASE ( 'single_street_canyon' )
448          ch  = NINT( canyon_height / dz )
449          IF ( canyon_width_x /= 9999999.9 )  THEN
450!
451!--          Street canyon in y direction
452             cwx = NINT( canyon_width_x / dx )
453             IF ( canyon_wall_left == 9999999.9 )  THEN
454                canyon_wall_left = ( nx + 1 - cwx ) / 2 * dx
455             ENDIF
456             cxl = NINT( canyon_wall_left / dx )
457             cxr = cxl + cwx
458             WRITE ( io, 272 )  'y', canyon_height, ch, 'u', cxl, cxr
459
460          ELSEIF ( canyon_width_y /= 9999999.9 )  THEN
461!
462!--          Street canyon in x direction
463             cwy = NINT( canyon_width_y / dy )
464             IF ( canyon_wall_south == 9999999.9 )  THEN
465                canyon_wall_south = ( ny + 1 - cwy ) / 2 * dy
466             ENDIF
467             cys = NINT( canyon_wall_south / dy )
468             cyn = cys + cwy
469             WRITE ( io, 272 )  'x', canyon_height, ch, 'v', cys, cyn
470          ENDIF
471
472    END SELECT
473
474    IF ( TRIM( topography ) /= 'flat' )  THEN
475       IF ( TRIM( topography_grid_convention ) == ' ' )  THEN
476          IF ( TRIM( topography ) == 'single_building' .OR.  &
477               TRIM( topography ) == 'single_street_canyon' )  THEN
478             WRITE ( io, 278 )
479          ELSEIF ( TRIM( topography ) == 'read_from_file' )  THEN
480             WRITE ( io, 279 )
481          ENDIF
482       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_edge' )  THEN
483          WRITE ( io, 278 )
484       ELSEIF ( TRIM( topography_grid_convention ) == 'cell_center' )  THEN
485          WRITE ( io, 279 )
486       ENDIF
487    ENDIF
488
489    IF ( plant_canopy ) THEN
490
491       WRITE ( io, 280 ) canopy_mode, pch_index, drag_coefficient
492       IF ( passive_scalar ) THEN
493          WRITE ( io, 281 ) scalar_exchange_coefficient,   &
494                            leaf_surface_concentration
495       ENDIF
496
497!
498!--    Heat flux at the top of vegetation
499       WRITE ( io, 282 ) cthf
500
501!
502!--    Leaf area density profile
503!--    Building output strings, starting with surface value
504       WRITE ( learde, '(F6.2)' )  lad_surface
505       gradients = '------'
506       slices = '     0'
507       coordinates = '   0.0'
508       i = 1
509       DO  WHILE ( lad_vertical_gradient_level_ind(i) /= -9999 )
510
511          WRITE (coor_chr,'(F7.2)')  lad(lad_vertical_gradient_level_ind(i))
512          learde = TRIM( learde ) // ' ' // TRIM( coor_chr )
513
514          WRITE (coor_chr,'(F7.2)')  lad_vertical_gradient(i)
515          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
516
517          WRITE (coor_chr,'(I7)')  lad_vertical_gradient_level_ind(i)
518          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
519
520          WRITE (coor_chr,'(F7.1)')  lad_vertical_gradient_level(i)
521          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
522
523          i = i + 1
524       ENDDO
525
526       WRITE ( io, 283 )  TRIM( coordinates ), TRIM( learde ), &
527                          TRIM( gradients ), TRIM( slices )
528
529    ENDIF
530
531!
532!-- Boundary conditions
533    IF ( ibc_p_b == 0 )  THEN
534       runten = 'p(0)     = 0      |'
535    ELSEIF ( ibc_p_b == 1 )  THEN
536       runten = 'p(0)     = p(1)   |'
537    ELSE
538       runten = 'p(0)     = p(1) +R|'
539    ENDIF
540    IF ( ibc_p_t == 0 )  THEN
541       roben  = 'p(nzt+1) = 0      |'
542    ELSE
543       roben  = 'p(nzt+1) = p(nzt) |'
544    ENDIF
545
546    IF ( ibc_uv_b == 0 )  THEN
547       runten = TRIM( runten ) // ' uv(0)     = -uv(1)                |'
548    ELSE
549       runten = TRIM( runten ) // ' uv(0)     = uv(1)                 |'
550    ENDIF
551    IF ( TRIM( bc_uv_t ) == 'dirichlet_0' )  THEN
552       roben  = TRIM( roben  ) // ' uv(nzt+1) = 0                     |'
553    ELSEIF ( ibc_uv_t == 0 )  THEN
554       roben  = TRIM( roben  ) // ' uv(nzt+1) = ug(nzt+1), vg(nzt+1)  |'
555    ELSE
556       roben  = TRIM( roben  ) // ' uv(nzt+1) = uv(nzt)               |'
557    ENDIF
558
559    IF ( ibc_pt_b == 0 )  THEN
560       runten = TRIM( runten ) // ' pt(0)   = pt_surface'
561    ELSEIF ( ibc_pt_b == 1 )  THEN
562       runten = TRIM( runten ) // ' pt(0)   = pt(1)'
563    ELSEIF ( ibc_pt_b == 2 )  THEN
564       runten = TRIM( runten ) // ' pt(0) = from coupled model'
565    ENDIF
566    IF ( ibc_pt_t == 0 )  THEN
567       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt_top'
568    ELSEIF( ibc_pt_t == 1 )  THEN
569       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt)'
570    ELSEIF( ibc_pt_t == 2 )  THEN
571       roben  = TRIM( roben  ) // ' pt(nzt+1) = pt(nzt) + dpt/dz_ini'
572    ENDIF
573
574    WRITE ( io, 300 )  runten, roben
575
576    IF ( .NOT. constant_diffusion )  THEN
577       IF ( ibc_e_b == 1 )  THEN
578          runten = 'e(0)     = e(1)'
579       ELSE
580          runten = 'e(0)     = e(1) = (u*/0.1)**2'
581       ENDIF
582       roben = 'e(nzt+1) = e(nzt) = e(nzt-1)'
583
584       WRITE ( io, 301 )  'e', runten, roben       
585
586    ENDIF
587
588    IF ( ocean )  THEN
589       runten = 'sa(0)    = sa(1)'
590       IF ( ibc_sa_t == 0 )  THEN
591          roben =  'sa(nzt+1) = sa_surface'
592       ELSE
593          roben =  'sa(nzt+1) = sa(nzt)'
594       ENDIF
595       WRITE ( io, 301 ) 'sa', runten, roben
596    ENDIF
597
598    IF ( humidity )  THEN
599       IF ( ibc_q_b == 0 )  THEN
600          runten = 'q(0)     = q_surface'
601       ELSE
602          runten = 'q(0)     = q(1)'
603       ENDIF
604       IF ( ibc_q_t == 0 )  THEN
605          roben =  'q(nzt)   = q_top'
606       ELSE
607          roben =  'q(nzt)   = q(nzt-1) + dq/dz'
608       ENDIF
609       WRITE ( io, 301 ) 'q', runten, roben
610    ENDIF
611
612    IF ( passive_scalar )  THEN
613       IF ( ibc_q_b == 0 )  THEN
614          runten = 's(0)     = s_surface'
615       ELSE
616          runten = 's(0)     = s(1)'
617       ENDIF
618       IF ( ibc_q_t == 0 )  THEN
619          roben =  's(nzt)   = s_top'
620       ELSE
621          roben =  's(nzt)   = s(nzt-1) + ds/dz'
622       ENDIF
623       WRITE ( io, 301 ) 's', runten, roben
624    ENDIF
625
626    IF ( use_surface_fluxes )  THEN
627       WRITE ( io, 303 )
628       IF ( constant_heatflux )  THEN
629          WRITE ( io, 306 )  surface_heatflux
630          IF ( random_heatflux )  WRITE ( io, 307 )
631       ENDIF
632       IF ( humidity  .AND.  constant_waterflux )  THEN
633          WRITE ( io, 311 ) surface_waterflux
634       ENDIF
635       IF ( passive_scalar  .AND.  constant_waterflux )  THEN
636          WRITE ( io, 313 ) surface_waterflux
637       ENDIF
638    ENDIF
639
640    IF ( use_top_fluxes )  THEN
641       WRITE ( io, 304 )
642       IF ( coupling_mode == 'uncoupled' )  THEN
643          WRITE ( io, 320 )  top_momentumflux_u, top_momentumflux_v
644          IF ( constant_top_heatflux )  THEN
645             WRITE ( io, 306 )  top_heatflux
646          ENDIF
647       ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
648          WRITE ( io, 316 )
649       ENDIF
650       IF ( ocean  .AND.  constant_top_salinityflux )  THEN
651          WRITE ( io, 309 )  top_salinityflux
652       ENDIF
653       IF ( humidity  .OR.  passive_scalar )  THEN
654          WRITE ( io, 315 )
655       ENDIF
656    ENDIF
657
658    IF ( prandtl_layer )  THEN
659       WRITE ( io, 305 )  0.5 * (zu(1)-zu(0)), roughness_length, kappa, &
660                          rif_min, rif_max
661       IF ( .NOT. constant_heatflux )  WRITE ( io, 308 )
662       IF ( humidity  .AND.  .NOT. constant_waterflux )  THEN
663          WRITE ( io, 312 )
664       ENDIF
665       IF ( passive_scalar  .AND.  .NOT. constant_waterflux )  THEN
666          WRITE ( io, 314 )
667       ENDIF
668    ELSE
669       IF ( INDEX(initializing_actions, 'set_1d-model_profiles') /= 0 )  THEN
670          WRITE ( io, 310 )  rif_min, rif_max
671       ENDIF
672    ENDIF
673
674    WRITE ( io, 317 )  bc_lr, bc_ns
675    IF ( bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic' )  THEN
676       WRITE ( io, 318 )  outflow_damping_width, km_damp_max
677       IF ( turbulent_inflow )  THEN
678          WRITE ( io, 319 )  recycling_width, recycling_plane, &
679                             inflow_damping_height, inflow_damping_width
680       ENDIF
681    ENDIF
682
683!
684!-- Listing of 1D-profiles
685    WRITE ( io, 325 )  dt_dopr_listing
686    IF ( averaging_interval_pr /= 0.0 )  THEN
687       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
688    ENDIF
689
690!
691!-- DATA output
692    WRITE ( io, 330 )
693    IF ( averaging_interval_pr /= 0.0 )  THEN
694       WRITE ( io, 326 )  averaging_interval_pr, dt_averaging_input_pr
695    ENDIF
696
697!
698!-- 1D-profiles
699    dopr_chr = 'Profile:'
700    IF ( dopr_n /= 0 )  THEN
701       WRITE ( io, 331 )
702
703       output_format = ''
704       IF ( netcdf_output )  THEN
705          IF ( netcdf_data_format == 1 )  THEN
706             output_format = 'NetCDF classic'
707          ELSE
708             output_format = 'NetCDF 64bit offset'
709          ENDIF
710       ENDIF
711       IF ( profil_output )  THEN
712          IF ( netcdf_output )  THEN
713             output_format = TRIM( output_format ) // ' and profil'
714          ELSE
715             output_format = 'profil'
716          ENDIF
717       ENDIF
718       WRITE ( io, 344 )  output_format
719
720       DO  i = 1, dopr_n
721          dopr_chr = TRIM( dopr_chr ) // ' ' // TRIM( data_output_pr(i) ) // ','
722          IF ( LEN_TRIM( dopr_chr ) >= 60 )  THEN
723             WRITE ( io, 332 )  dopr_chr
724             dopr_chr = '       :'
725          ENDIF
726       ENDDO
727
728       IF ( dopr_chr /= '' )  THEN
729          WRITE ( io, 332 )  dopr_chr
730       ENDIF
731       WRITE ( io, 333 )  dt_dopr, averaging_interval_pr, dt_averaging_input_pr
732       IF ( skip_time_dopr /= 0.0 )  WRITE ( io, 339 )  skip_time_dopr
733    ENDIF
734
735!
736!-- 2D-arrays
737    DO  av = 0, 1
738
739       i = 1
740       do2d_xy = ''
741       do2d_xz = ''
742       do2d_yz = ''
743       DO  WHILE ( do2d(av,i) /= ' ' )
744
745          l = MAX( 2, LEN_TRIM( do2d(av,i) ) )
746          do2d_mode = do2d(av,i)(l-1:l)
747
748          SELECT CASE ( do2d_mode )
749             CASE ( 'xy' )
750                ll = LEN_TRIM( do2d_xy )
751                do2d_xy = do2d_xy(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
752             CASE ( 'xz' )
753                ll = LEN_TRIM( do2d_xz )
754                do2d_xz = do2d_xz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
755             CASE ( 'yz' )
756                ll = LEN_TRIM( do2d_yz )
757                do2d_yz = do2d_yz(1:ll) // ' ' // do2d(av,i)(1:l-3) // ','
758          END SELECT
759
760          i = i + 1
761
762       ENDDO
763
764       IF ( ( ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  .OR.    &
765              ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  .OR.    &
766              ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 ) )  .AND. &
767            ( netcdf_output  .OR.  iso2d_output ) )  THEN
768
769          IF (  av == 0 )  THEN
770             WRITE ( io, 334 )  ''
771          ELSE
772             WRITE ( io, 334 )  '(time-averaged)'
773          ENDIF
774
775          IF ( do2d_at_begin )  THEN
776             begin_chr = 'and at the start'
777          ELSE
778             begin_chr = ''
779          ENDIF
780
781          output_format = ''
782          IF ( netcdf_output )  THEN
783             IF ( netcdf_data_format == 1 )  THEN
784                output_format = 'NetCDF classic'
785             ELSEIF ( netcdf_data_format == 2 )  THEN
786                output_format = 'NetCDF 64bit offset'
787             ELSEIF ( netcdf_data_format == 3 )  THEN
788                output_format = 'NetCDF4/HDF5'
789             ELSEIF ( netcdf_data_format == 4 )  THEN
790                output_format = 'NetCDF4/HDF5 clasic'
791             ENDIF
792          ENDIF
793          IF ( iso2d_output )  THEN
794             IF ( netcdf_output )  THEN
795                output_format = TRIM( output_format ) // ' and iso2d'
796             ELSE
797                output_format = 'iso2d'
798             ENDIF
799          ENDIF
800          WRITE ( io, 344 )  output_format
801
802          IF ( do2d_xy /= ''  .AND.  section(1,1) /= -9999 )  THEN
803             i = 1
804             slices = '/'
805             coordinates = '/'
806!
807!--          Building strings with index and coordinate informations of the
808!--          slices
809             DO  WHILE ( section(i,1) /= -9999 )
810
811                WRITE (section_chr,'(I5)')  section(i,1)
812                section_chr = ADJUSTL( section_chr )
813                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
814
815                IF ( section(i,1) == -1 )  THEN
816                   WRITE (coor_chr,'(F10.1)')  -1.0
817                ELSE
818                   WRITE (coor_chr,'(F10.1)')  zu(section(i,1))
819                ENDIF
820                coor_chr = ADJUSTL( coor_chr )
821                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
822
823                i = i + 1
824             ENDDO
825             IF ( av == 0 )  THEN
826                WRITE ( io, 335 )  'XY', do2d_xy, dt_do2d_xy, &
827                                   TRIM( begin_chr ), 'k', TRIM( slices ), &
828                                   TRIM( coordinates )
829                IF ( skip_time_do2d_xy /= 0.0 )  THEN
830                   WRITE ( io, 339 )  skip_time_do2d_xy
831                ENDIF
832             ELSE
833                WRITE ( io, 342 )  'XY', do2d_xy, dt_data_output_av, &
834                                   TRIM( begin_chr ), averaging_interval, &
835                                   dt_averaging_input, 'k', TRIM( slices ), &
836                                   TRIM( coordinates )
837                IF ( skip_time_data_output_av /= 0.0 )  THEN
838                   WRITE ( io, 339 )  skip_time_data_output_av
839                ENDIF
840             ENDIF
841
842          ENDIF
843
844          IF ( do2d_xz /= ''  .AND.  section(1,2) /= -9999 )  THEN
845             i = 1
846             slices = '/'
847             coordinates = '/'
848!
849!--          Building strings with index and coordinate informations of the
850!--          slices
851             DO  WHILE ( section(i,2) /= -9999 )
852
853                WRITE (section_chr,'(I5)')  section(i,2)
854                section_chr = ADJUSTL( section_chr )
855                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
856
857                WRITE (coor_chr,'(F10.1)')  section(i,2) * dy
858                coor_chr = ADJUSTL( coor_chr )
859                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
860
861                i = i + 1
862             ENDDO
863             IF ( av == 0 )  THEN
864                WRITE ( io, 335 )  'XZ', do2d_xz, dt_do2d_xz, &
865                                   TRIM( begin_chr ), 'j', TRIM( slices ), &
866                                   TRIM( coordinates )
867                IF ( skip_time_do2d_xz /= 0.0 )  THEN
868                   WRITE ( io, 339 )  skip_time_do2d_xz
869                ENDIF
870             ELSE
871                WRITE ( io, 342 )  'XZ', do2d_xz, dt_data_output_av, &
872                                   TRIM( begin_chr ), averaging_interval, &
873                                   dt_averaging_input, 'j', TRIM( slices ), &
874                                   TRIM( coordinates )
875                IF ( skip_time_data_output_av /= 0.0 )  THEN
876                   WRITE ( io, 339 )  skip_time_data_output_av
877                ENDIF
878             ENDIF
879          ENDIF
880
881          IF ( do2d_yz /= ''  .AND.  section(1,3) /= -9999 )  THEN
882             i = 1
883             slices = '/'
884             coordinates = '/'
885!
886!--          Building strings with index and coordinate informations of the
887!--          slices
888             DO  WHILE ( section(i,3) /= -9999 )
889
890                WRITE (section_chr,'(I5)')  section(i,3)
891                section_chr = ADJUSTL( section_chr )
892                slices = TRIM( slices ) // TRIM( section_chr ) // '/'
893
894                WRITE (coor_chr,'(F10.1)')  section(i,3) * dx
895                coor_chr = ADJUSTL( coor_chr )
896                coordinates = TRIM( coordinates ) // TRIM( coor_chr ) // '/'
897
898                i = i + 1
899             ENDDO
900             IF ( av == 0 )  THEN
901                WRITE ( io, 335 )  'YZ', do2d_yz, dt_do2d_yz, &
902                                   TRIM( begin_chr ), 'i', TRIM( slices ), &
903                                   TRIM( coordinates )
904                IF ( skip_time_do2d_yz /= 0.0 )  THEN
905                   WRITE ( io, 339 )  skip_time_do2d_yz
906                ENDIF
907             ELSE
908                WRITE ( io, 342 )  'YZ', do2d_yz, dt_data_output_av, &
909                                   TRIM( begin_chr ), averaging_interval, &
910                                   dt_averaging_input, 'i', TRIM( slices ), &
911                                   TRIM( coordinates )
912                IF ( skip_time_data_output_av /= 0.0 )  THEN
913                   WRITE ( io, 339 )  skip_time_data_output_av
914                ENDIF
915             ENDIF
916          ENDIF
917
918       ENDIF
919
920    ENDDO
921
922!
923!-- 3d-arrays
924    DO  av = 0, 1
925
926       i = 1
927       do3d_chr = ''
928       DO  WHILE ( do3d(av,i) /= ' ' )
929
930          do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // ','
931          i = i + 1
932
933       ENDDO
934
935       IF ( do3d_chr /= '' )  THEN
936          IF ( av == 0 )  THEN
937             WRITE ( io, 336 )  ''
938          ELSE
939             WRITE ( io, 336 )  '(time-averaged)'
940          ENDIF
941
942          output_format = ''
943          IF ( netcdf_output )  THEN
944             IF ( netcdf_data_format == 1 )  THEN
945                output_format = 'NetCDF classic'
946             ELSEIF ( netcdf_data_format == 2 )  THEN
947                output_format = 'NetCDF 64bit offset'
948             ELSEIF ( netcdf_data_format == 3 )  THEN
949                output_format = 'NetCDF4/HDF5'
950             ELSEIF ( netcdf_data_format == 4 )  THEN
951                output_format = 'NetCDF4/HDF5 clasic'
952             ENDIF
953          ENDIF
954          IF ( avs_output )  THEN
955             IF ( netcdf_output )  THEN
956                output_format = TRIM( output_format ) // ' and avs'
957             ELSE
958                output_format = 'avs'
959             ENDIF
960          ENDIF
961          WRITE ( io, 344 )  output_format
962
963          IF ( do3d_at_begin )  THEN
964             begin_chr = 'and at the start'
965          ELSE
966             begin_chr = ''
967          ENDIF
968          IF ( av == 0 )  THEN
969             WRITE ( io, 337 )  do3d_chr, dt_do3d, TRIM( begin_chr ), &
970                                zu(nz_do3d), nz_do3d
971          ELSE
972             WRITE ( io, 343 )  do3d_chr, dt_data_output_av,           &
973                                TRIM( begin_chr ), averaging_interval, &
974                                dt_averaging_input, zu(nz_do3d), nz_do3d
975          ENDIF
976
977          IF ( do3d_compress )  THEN
978             do3d_chr = ''
979             i = 1
980             DO WHILE ( do3d(av,i) /= ' ' )
981
982                SELECT CASE ( do3d(av,i) )
983                   CASE ( 'u' )
984                      j = 1
985                   CASE ( 'v' )
986                      j = 2
987                   CASE ( 'w' )
988                      j = 3
989                   CASE ( 'p' )
990                      j = 4
991                   CASE ( 'pt' )
992                      j = 5
993                END SELECT
994                WRITE ( prec, '(I1)' )  plot_3d_precision(j)%precision
995                do3d_chr = TRIM( do3d_chr ) // ' ' // TRIM( do3d(av,i) ) // &
996                           ':' // prec // ','
997                i = i + 1
998
999             ENDDO
1000             WRITE ( io, 338 )  do3d_chr
1001
1002          ENDIF
1003
1004          IF ( av == 0 )  THEN
1005             IF ( skip_time_do3d /= 0.0 )  THEN
1006                WRITE ( io, 339 )  skip_time_do3d
1007             ENDIF
1008          ELSE
1009             IF ( skip_time_data_output_av /= 0.0 )  THEN
1010                WRITE ( io, 339 )  skip_time_data_output_av
1011             ENDIF
1012          ENDIF
1013
1014       ENDIF
1015
1016    ENDDO
1017
1018!
1019!-- masked arrays
1020    IF ( masks > 0 )  WRITE ( io, 345 )  &
1021         mask_scale_x, mask_scale_y, mask_scale_z
1022    DO  mid = 1, masks
1023       DO  av = 0, 1
1024
1025          i = 1
1026          domask_chr = ''
1027          DO  WHILE ( domask(mid,av,i) /= ' ' )
1028             domask_chr = TRIM( domask_chr ) // ' ' //  &
1029                          TRIM( domask(mid,av,i) ) // ','
1030             i = i + 1
1031          ENDDO
1032
1033          IF ( domask_chr /= '' )  THEN
1034             IF ( av == 0 )  THEN
1035                WRITE ( io, 346 )  '', mid
1036             ELSE
1037                WRITE ( io, 346 )  ' (time-averaged)', mid
1038             ENDIF
1039
1040             output_format = ''
1041             IF ( netcdf_output )  THEN
1042                IF ( netcdf_data_format == 1 )  THEN
1043                   output_format = 'NetCDF classic'
1044                ELSEIF ( netcdf_data_format == 2 )  THEN
1045                   output_format = 'NetCDF 64bit offset'
1046                ELSEIF ( netcdf_data_format == 3 )  THEN
1047                   output_format = 'NetCDF4/HDF5'
1048                ELSEIF ( netcdf_data_format == 4 )  THEN
1049                   output_format = 'NetCDF4/HDF5 clasic'
1050                ENDIF
1051             ENDIF
1052             WRITE ( io, 344 )  output_format
1053
1054             IF ( av == 0 )  THEN
1055                WRITE ( io, 347 )  domask_chr, dt_domask(mid)
1056             ELSE
1057                WRITE ( io, 348 )  domask_chr, dt_data_output_av, &
1058                                   averaging_interval, dt_averaging_input
1059             ENDIF
1060
1061             IF ( av == 0 )  THEN
1062                IF ( skip_time_domask(mid) /= 0.0 )  THEN
1063                   WRITE ( io, 339 )  skip_time_domask(mid)
1064                ENDIF
1065             ELSE
1066                IF ( skip_time_data_output_av /= 0.0 )  THEN
1067                   WRITE ( io, 339 )  skip_time_data_output_av
1068                ENDIF
1069             ENDIF
1070!
1071!--          output locations
1072             DO  dim = 1, 3
1073                IF ( mask(mid,dim,1) >= 0.0 )  THEN
1074                   count = 0
1075                   DO  WHILE ( mask(mid,dim,count+1) >= 0.0 )
1076                      count = count + 1
1077                   ENDDO
1078                   WRITE ( io, 349 )  dir(dim), dir(dim), mid, dir(dim), &
1079                                      mask(mid,dim,:count)
1080                ELSEIF ( mask_loop(mid,dim,1) < 0.0 .AND.  &
1081                         mask_loop(mid,dim,2) < 0.0 .AND.  &
1082                         mask_loop(mid,dim,3) == 0.0 )  THEN
1083                   WRITE ( io, 350 )  dir(dim), dir(dim)
1084                ELSEIF ( mask_loop(mid,dim,3) == 0.0 )  THEN
1085                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1086                                      mask_loop(mid,dim,1:2)
1087                ELSE
1088                   WRITE ( io, 351 )  dir(dim), dir(dim), mid, dir(dim), &
1089                                      mask_loop(mid,dim,1:3)
1090                ENDIF
1091             ENDDO
1092          ENDIF
1093
1094       ENDDO
1095    ENDDO
1096
1097!
1098!-- Timeseries
1099    IF ( dt_dots /= 9999999.9 )  THEN
1100       WRITE ( io, 340 )
1101
1102       output_format = ''
1103       IF ( netcdf_output )  THEN
1104          IF ( netcdf_data_format == 1 )  THEN
1105             output_format = 'NetCDF classic'
1106          ELSE
1107             output_format = 'NetCDF 64bit offset'
1108          ENDIF
1109       ENDIF
1110       IF ( profil_output )  THEN
1111          IF ( netcdf_output )  THEN
1112             output_format = TRIM( output_format ) // ' and profil'
1113          ELSE
1114             output_format = 'profil'
1115          ENDIF
1116       ENDIF
1117       WRITE ( io, 344 )  output_format
1118       WRITE ( io, 341 )  dt_dots
1119    ENDIF
1120
1121#if defined( __dvrp_graphics )
1122!
1123!-- Dvrp-output
1124    IF ( dt_dvrp /= 9999999.9 )  THEN
1125       WRITE ( io, 360 )  dt_dvrp, TRIM( dvrp_output ), TRIM( dvrp_host ), &
1126                          TRIM( dvrp_username ), TRIM( dvrp_directory )
1127       i = 1
1128       l = 0
1129       m = 0
1130       DO WHILE ( mode_dvrp(i) /= ' ' )
1131          IF ( mode_dvrp(i)(1:10) == 'isosurface' )  THEN
1132             READ ( mode_dvrp(i), '(10X,I2)' )  j
1133             l = l + 1
1134             IF ( do3d(0,j) /= ' ' )  THEN
1135                WRITE ( io, 361 )  TRIM( do3d(0,j) ), threshold(l), &
1136                                   isosurface_color(:,l)
1137             ENDIF
1138          ELSEIF ( mode_dvrp(i)(1:6) == 'slicer' )  THEN
1139             READ ( mode_dvrp(i), '(6X,I2)' )  j
1140             m = m + 1
1141             IF ( do2d(0,j) /= ' ' )  THEN
1142                WRITE ( io, 362 )  TRIM( do2d(0,j) ), &
1143                                   slicer_range_limits_dvrp(:,m)
1144             ENDIF
1145          ELSEIF ( mode_dvrp(i)(1:9) == 'particles' )  THEN
1146             WRITE ( io, 363 )  dvrp_psize
1147             IF ( particle_dvrpsize /= 'none' )  THEN
1148                WRITE ( io, 364 )  'size', TRIM( particle_dvrpsize ), &
1149                                   dvrpsize_interval
1150             ENDIF
1151             IF ( particle_color /= 'none' )  THEN
1152                WRITE ( io, 364 )  'color', TRIM( particle_color ), &
1153                                   color_interval
1154             ENDIF
1155          ENDIF
1156          i = i + 1
1157       ENDDO
1158
1159       WRITE ( io, 365 )  groundplate_color, superelevation_x, &
1160                          superelevation_y, superelevation, clip_dvrp_l, &
1161                          clip_dvrp_r, clip_dvrp_s, clip_dvrp_n
1162
1163       IF ( TRIM( topography ) /= 'flat' )  THEN
1164          WRITE ( io, 366 )  topography_color
1165          IF ( cluster_size > 1 )  THEN
1166             WRITE ( io, 367 )  cluster_size
1167          ENDIF
1168       ENDIF
1169
1170    ENDIF
1171#endif
1172
1173#if defined( __spectra )
1174!
1175!-- Spectra output
1176    IF ( dt_dosp /= 9999999.9 ) THEN
1177       WRITE ( io, 370 )
1178
1179       output_format = ''
1180       IF ( netcdf_output )  THEN
1181          IF ( netcdf_data_format == 1 )  THEN
1182             output_format = 'NetCDF classic'
1183          ELSE
1184             output_format = 'NetCDF 64bit offset'
1185          ENDIF
1186       ENDIF
1187       IF ( profil_output )  THEN
1188          IF ( netcdf_output )  THEN
1189             output_format = TRIM( output_format ) // ' and profil'
1190          ELSE
1191             output_format = 'profil'
1192          ENDIF
1193       ENDIF
1194       WRITE ( io, 344 )  output_format
1195       WRITE ( io, 371 )  dt_dosp
1196       IF ( skip_time_dosp /= 0.0 )  WRITE ( io, 339 )  skip_time_dosp
1197       WRITE ( io, 372 )  ( data_output_sp(i), i = 1,10 ),     &
1198                          ( spectra_direction(i), i = 1,10 ),  &
1199                          ( comp_spectra_level(i), i = 1,100 ), &
1200                          ( plot_spectra_level(i), i = 1,100 ), &
1201                          averaging_interval_sp, dt_averaging_input_pr
1202    ENDIF
1203#endif
1204
1205    WRITE ( io, 99 )
1206
1207!
1208!-- Physical quantities
1209    WRITE ( io, 400 )
1210
1211!
1212!-- Geostrophic parameters
1213    WRITE ( io, 410 )  omega, phi, f, fs
1214
1215!
1216!-- Other quantities
1217    WRITE ( io, 411 )  g
1218    IF ( use_reference )  THEN
1219       IF ( ocean )  THEN
1220          WRITE ( io, 412 )  prho_reference
1221       ELSE
1222          WRITE ( io, 413 )  pt_reference
1223       ENDIF
1224    ENDIF
1225
1226!
1227!-- Cloud physics parameters
1228    IF ( cloud_physics ) THEN
1229       WRITE ( io, 415 )
1230       WRITE ( io, 416 ) surface_pressure, r_d, rho_surface, cp, l_v
1231    ENDIF
1232
1233!-- Profile of the geostrophic wind (component ug)
1234!-- Building output strings
1235    WRITE ( ugcomponent, '(F6.2)' )  ug_surface
1236    gradients = '------'
1237    slices = '     0'
1238    coordinates = '   0.0'
1239    i = 1
1240    DO  WHILE ( ug_vertical_gradient_level_ind(i) /= -9999 )
1241     
1242       WRITE (coor_chr,'(F6.2,1X)')  ug(ug_vertical_gradient_level_ind(i))
1243       ugcomponent = TRIM( ugcomponent ) // '  ' // TRIM( coor_chr )
1244
1245       WRITE (coor_chr,'(F6.2,1X)')  ug_vertical_gradient(i)
1246       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1247
1248       WRITE (coor_chr,'(I6,1X)')  ug_vertical_gradient_level_ind(i)
1249       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1250
1251       WRITE (coor_chr,'(F6.1,1X)')  ug_vertical_gradient_level(i)
1252       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
1253
1254       IF ( i == 10 )  THEN
1255          EXIT
1256       ELSE
1257          i = i + 1
1258       ENDIF
1259
1260    ENDDO
1261
1262    WRITE ( io, 423 )  TRIM( coordinates ), TRIM( ugcomponent ), &
1263                       TRIM( gradients ), TRIM( slices )
1264
1265!-- Profile of the geostrophic wind (component vg)
1266!-- Building output strings
1267    WRITE ( vgcomponent, '(F6.2)' )  vg_surface
1268    gradients = '------'
1269    slices = '     0'
1270    coordinates = '   0.0'
1271    i = 1
1272    DO  WHILE ( vg_vertical_gradient_level_ind(i) /= -9999 )
1273
1274       WRITE (coor_chr,'(F6.2,1X)')  vg(vg_vertical_gradient_level_ind(i))
1275       vgcomponent = TRIM( vgcomponent ) // '  ' // TRIM( coor_chr )
1276
1277       WRITE (coor_chr,'(F6.2,1X)')  vg_vertical_gradient(i)
1278       gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1279
1280       WRITE (coor_chr,'(I6,1X)')  vg_vertical_gradient_level_ind(i)
1281       slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1282
1283       WRITE (coor_chr,'(F6.1,1X)')  vg_vertical_gradient_level(i)
1284       coordinates = TRIM( coordinates ) // '  ' // TRIM( coor_chr )
1285
1286       IF ( i == 10 )  THEN
1287          EXIT
1288       ELSE
1289          i = i + 1
1290       ENDIF
1291 
1292    ENDDO
1293
1294    WRITE ( io, 424 )  TRIM( coordinates ), TRIM( vgcomponent ), &
1295                       TRIM( gradients ), TRIM( slices )
1296
1297!
1298!-- Initial temperature profile
1299!-- Building output strings, starting with surface temperature
1300    WRITE ( temperatures, '(F6.2)' )  pt_surface
1301    gradients = '------'
1302    slices = '     0'
1303    coordinates = '   0.0'
1304    i = 1
1305    DO  WHILE ( pt_vertical_gradient_level_ind(i) /= -9999 )
1306
1307       WRITE (coor_chr,'(F7.2)')  pt_init(pt_vertical_gradient_level_ind(i))
1308       temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1309
1310       WRITE (coor_chr,'(F7.2)')  pt_vertical_gradient(i)
1311       gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1312
1313       WRITE (coor_chr,'(I7)')  pt_vertical_gradient_level_ind(i)
1314       slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1315
1316       WRITE (coor_chr,'(F7.1)')  pt_vertical_gradient_level(i)
1317       coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1318
1319       IF ( i == 10 )  THEN
1320          EXIT
1321       ELSE
1322          i = i + 1
1323       ENDIF
1324
1325    ENDDO
1326
1327    WRITE ( io, 420 )  TRIM( coordinates ), TRIM( temperatures ), &
1328                       TRIM( gradients ), TRIM( slices )
1329
1330!
1331!-- Initial humidity profile
1332!-- Building output strings, starting with surface humidity
1333    IF ( humidity  .OR.  passive_scalar )  THEN
1334       WRITE ( temperatures, '(E8.1)' )  q_surface
1335       gradients = '--------'
1336       slices = '       0'
1337       coordinates = '     0.0'
1338       i = 1
1339       DO  WHILE ( q_vertical_gradient_level_ind(i) /= -9999 )
1340         
1341          WRITE (coor_chr,'(E8.1,4X)')  q_init(q_vertical_gradient_level_ind(i))
1342          temperatures = TRIM( temperatures ) // '  ' // TRIM( coor_chr )
1343
1344          WRITE (coor_chr,'(E8.1,4X)')  q_vertical_gradient(i)
1345          gradients = TRIM( gradients ) // '  ' // TRIM( coor_chr )
1346         
1347          WRITE (coor_chr,'(I8,4X)')  q_vertical_gradient_level_ind(i)
1348          slices = TRIM( slices ) // '  ' // TRIM( coor_chr )
1349         
1350          WRITE (coor_chr,'(F8.1,4X)')  q_vertical_gradient_level(i)
1351          coordinates = TRIM( coordinates ) // '  '  // TRIM( coor_chr )
1352
1353          IF ( i == 10 )  THEN
1354             EXIT
1355          ELSE
1356             i = i + 1
1357          ENDIF
1358
1359       ENDDO
1360
1361       IF ( humidity )  THEN
1362          WRITE ( io, 421 )  TRIM( coordinates ), TRIM( temperatures ), &
1363                             TRIM( gradients ), TRIM( slices )
1364       ELSE
1365          WRITE ( io, 422 )  TRIM( coordinates ), TRIM( temperatures ), &
1366                             TRIM( gradients ), TRIM( slices )
1367       ENDIF
1368    ENDIF
1369
1370!
1371!-- Initial salinity profile
1372!-- Building output strings, starting with surface salinity
1373    IF ( ocean )  THEN
1374       WRITE ( temperatures, '(F6.2)' )  sa_surface
1375       gradients = '------'
1376       slices = '     0'
1377       coordinates = '   0.0'
1378       i = 1
1379       DO  WHILE ( sa_vertical_gradient_level_ind(i) /= -9999 )
1380
1381          WRITE (coor_chr,'(F7.2)')  sa_init(sa_vertical_gradient_level_ind(i))
1382          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1383
1384          WRITE (coor_chr,'(F7.2)')  sa_vertical_gradient(i)
1385          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1386
1387          WRITE (coor_chr,'(I7)')  sa_vertical_gradient_level_ind(i)
1388          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1389
1390          WRITE (coor_chr,'(F7.1)')  sa_vertical_gradient_level(i)
1391          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1392
1393          IF ( i == 10 )  THEN
1394             EXIT
1395          ELSE
1396             i = i + 1
1397          ENDIF
1398
1399       ENDDO
1400
1401       WRITE ( io, 425 )  TRIM( coordinates ), TRIM( temperatures ), &
1402                          TRIM( gradients ), TRIM( slices )
1403    ENDIF
1404
1405!
1406!-- Profile for the large scale vertial velocity
1407!-- Building output strings, starting with surface value
1408    IF ( large_scale_subsidence )  THEN
1409       temperatures = '   0.0'
1410       gradients = '------'
1411       slices = '     0'
1412       coordinates = '   0.0'
1413       i = 1
1414       DO  WHILE ( ws_vertical_gradient_level_ind(i) /= -9999 )
1415
1416          WRITE (coor_chr,'(E10.2,7X)')  &
1417                                w_subs(ws_vertical_gradient_level_ind(i))
1418          temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr )
1419
1420          WRITE (coor_chr,'(E10.2,7X)')  ws_vertical_gradient(i)
1421          gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr )
1422
1423          WRITE (coor_chr,'(I10,7X)')  ws_vertical_gradient_level_ind(i)
1424          slices = TRIM( slices ) // ' ' // TRIM( coor_chr )
1425
1426          WRITE (coor_chr,'(F10.2,7X)')  ws_vertical_gradient_level(i)
1427          coordinates = TRIM( coordinates ) // ' '  // TRIM( coor_chr )
1428
1429          IF ( i == 10 )  THEN
1430             EXIT
1431          ELSE
1432             i = i + 1
1433          ENDIF
1434
1435       ENDDO
1436
1437       WRITE ( io, 426 )  TRIM( coordinates ), TRIM( temperatures ), &
1438                          TRIM( gradients ), TRIM( slices )
1439    ENDIF
1440
1441!
1442!-- LES / turbulence parameters
1443    WRITE ( io, 450 )
1444
1445!--
1446! ... LES-constants used must still be added here
1447!--
1448    IF ( constant_diffusion )  THEN
1449       WRITE ( io, 451 )  km_constant, km_constant/prandtl_number, &
1450                          prandtl_number
1451    ENDIF
1452    IF ( .NOT. constant_diffusion)  THEN
1453       IF ( e_init > 0.0 )  WRITE ( io, 455 )  e_init
1454       IF ( e_min > 0.0 )  WRITE ( io, 454 )  e_min
1455       IF ( wall_adjustment )  WRITE ( io, 453 )  wall_adjustment_factor
1456       IF ( adjust_mixing_length  .AND.  prandtl_layer )  WRITE ( io, 452 )
1457    ENDIF
1458
1459!
1460!-- Special actions during the run
1461    WRITE ( io, 470 )
1462    IF ( create_disturbances )  THEN
1463       WRITE ( io, 471 )  dt_disturb, disturbance_amplitude,                   &
1464                          zu(disturbance_level_ind_b), disturbance_level_ind_b,&
1465                          zu(disturbance_level_ind_t), disturbance_level_ind_t
1466       IF ( bc_lr /= 'cyclic'  .OR.  bc_ns /= 'cyclic' )  THEN
1467          WRITE ( io, 472 )  inflow_disturbance_begin, inflow_disturbance_end
1468       ELSE
1469          WRITE ( io, 473 )  disturbance_energy_limit
1470       ENDIF
1471       WRITE ( io, 474 )  TRIM( random_generator )
1472    ENDIF
1473    IF ( pt_surface_initial_change /= 0.0 )  THEN
1474       WRITE ( io, 475 )  pt_surface_initial_change
1475    ENDIF
1476    IF ( humidity  .AND.  q_surface_initial_change /= 0.0 )  THEN
1477       WRITE ( io, 476 )  q_surface_initial_change       
1478    ENDIF
1479    IF ( passive_scalar  .AND.  q_surface_initial_change /= 0.0 )  THEN
1480       WRITE ( io, 477 )  q_surface_initial_change       
1481    ENDIF
1482
1483    IF ( particle_advection )  THEN
1484!
1485!--    Particle attributes
1486       WRITE ( io, 480 )  particle_advection_start, dt_prel, bc_par_lr, &
1487                          bc_par_ns, bc_par_b, bc_par_t, particle_maximum_age, &
1488                          end_time_prel, dt_sort_particles
1489       IF ( use_sgs_for_particles )  WRITE ( io, 488 )  dt_min_part
1490       IF ( random_start_position )  WRITE ( io, 481 )
1491       IF ( particles_per_point > 1 )  WRITE ( io, 489 )  particles_per_point
1492       WRITE ( io, 495 )  total_number_of_particles
1493       IF ( maximum_number_of_tailpoints /= 0 )  THEN
1494          WRITE ( io, 483 )  maximum_number_of_tailpoints
1495          IF ( minimum_tailpoint_distance /= 0 )  THEN
1496             WRITE ( io, 484 )  total_number_of_tails,      &
1497                                minimum_tailpoint_distance, &
1498                                maximum_tailpoint_age
1499          ENDIF
1500       ENDIF
1501       IF ( dt_write_particle_data /= 9999999.9 )  THEN
1502          WRITE ( io, 485 )  dt_write_particle_data
1503          output_format = ''
1504          IF ( netcdf_output )  THEN
1505             IF ( netcdf_data_format > 1 )  THEN
1506                output_format = 'netcdf (64 bit offset) and binary'
1507             ELSE
1508                output_format = 'netcdf and binary'
1509             ENDIF
1510          ELSE
1511             output_format = 'binary'
1512          ENDIF
1513          WRITE ( io, 344 )  output_format
1514       ENDIF
1515       IF ( dt_dopts /= 9999999.9 )  WRITE ( io, 494 )  dt_dopts
1516       IF ( write_particle_statistics )  WRITE ( io, 486 )
1517
1518       WRITE ( io, 487 )  number_of_particle_groups
1519
1520       DO  i = 1, number_of_particle_groups
1521          IF ( i == 1  .AND.  density_ratio(i) == 9999999.9 )  THEN
1522             WRITE ( io, 490 )  i, 0.0
1523             WRITE ( io, 492 )
1524          ELSE
1525             WRITE ( io, 490 )  i, radius(i)
1526             IF ( density_ratio(i) /= 0.0 )  THEN
1527                WRITE ( io, 491 )  density_ratio(i)
1528             ELSE
1529                WRITE ( io, 492 )
1530             ENDIF
1531          ENDIF
1532          WRITE ( io, 493 )  psl(i), psr(i), pss(i), psn(i), psb(i), pst(i), &
1533                             pdx(i), pdy(i), pdz(i)
1534          IF ( .NOT. vertical_particle_advection(i) )  WRITE ( io, 482 )
1535       ENDDO
1536
1537    ENDIF
1538
1539
1540!
1541!-- Parameters of 1D-model
1542    IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )  THEN
1543       WRITE ( io, 500 )  end_time_1d, dt_run_control_1d, dt_pr_1d, &
1544                          mixing_length_1d, dissipation_1d
1545       IF ( damp_level_ind_1d /= nzt+1 )  THEN
1546          WRITE ( io, 502 )  zu(damp_level_ind_1d), damp_level_ind_1d
1547       ENDIF
1548    ENDIF
1549
1550!
1551!-- User-defined informations
1552    CALL user_header( io )
1553
1554    WRITE ( io, 99 )
1555
1556!
1557!-- Write buffer contents to disc immediately
1558    CALL local_flush( io )
1559
1560!
1561!-- Here the FORMATs start
1562
1563 99 FORMAT (1X,78('-'))
1564100 FORMAT (/1X,'***************************',9X,42('-')/        &
1565            1X,'* ',A,' *',9X,A/                               &
1566            1X,'***************************',9X,42('-'))
1567101 FORMAT (37X,'coupled run using MPI-',I1,': ',A/ &
1568            37X,42('-'))
1569102 FORMAT (/' Date:              ',A8,9X,'Run:       ',A20/      &
1570            ' Time:              ',A8,9X,'Run-No.:   ',I2.2/     &
1571            ' Run on host:     ',A10)
1572#if defined( __parallel )
1573103 FORMAT (' Number of PEs:',8X,I5,9X,'Processor grid (x,y): (',I3,',',I3, &
1574              ')',1X,A)
1575104 FORMAT (' Number of PEs:',8X,I5,9X,'Tasks:',I4,'   threads per task:',I4/ &
1576              37X,'Processor grid (x,y): (',I3,',',I3,')',1X,A)
1577105 FORMAT (37X,'One additional PE is used to handle'/37X,'the dvrp output!')
1578106 FORMAT (37X,'A 1d-decomposition along x is forced'/ &
1579            37X,'because the job is running on an SMP-cluster')
1580107 FORMAT (37X,'A 1d-decomposition along ',A,' is used')
1581#endif
1582110 FORMAT (/' Numerical Schemes:'/ &
1583             ' -----------------'/)
1584111 FORMAT (' --> Solve perturbation pressure via FFT using ',A,' routines')
1585112 FORMAT (' --> Solve perturbation pressure via SOR-Red/Black-Schema'/ &
1586            '     Iterations (initial/other): ',I3,'/',I3,'  omega = ',F5.3)
1587113 FORMAT (' --> Momentum advection via Piascek-Williams-Scheme (Form C3)', &
1588                  ' or Upstream')
1589114 FORMAT (' --> Momentum advection via Upstream-Spline-Scheme')
1590115 FORMAT ('     Tendencies are smoothed via Long-Filter with factor ',F5.3) 
1591116 FORMAT (' --> Scalar advection via Piascek-Williams-Scheme (Form C3)', &
1592                  ' or Upstream')
1593117 FORMAT (' --> Scalar advection via Upstream-Spline-Scheme')
1594118 FORMAT (' --> Scalar advection via Bott-Chlond-Scheme')
1595119 FORMAT (' --> Galilei-Transform applied to horizontal advection', &
1596            '     Translation velocity = ',A/ &
1597            '     distance advected ',A,':  ',F8.3,' km(x)  ',F8.3,' km(y)')
1598120 FORMAT (' --> Time differencing scheme: leapfrog only (no euler in case', &
1599                  ' of timestep changes)')
1600121 FORMAT (' --> Time differencing scheme: leapfrog + euler in case of', &
1601                  ' timestep changes')
1602122 FORMAT (' --> Time differencing scheme: ',A)
1603123 FORMAT (' --> Rayleigh-Damping active, starts ',A,' z = ',F8.2,' m'/ &
1604            '     maximum damping coefficient: ',F5.3, ' 1/s')
1605124 FORMAT ('     Spline-overshoots are being suppressed')
1606125 FORMAT ('     Upstream-Scheme is used if Upstream-differences fall short', &
1607                  ' of'/                                                       &
1608            '     delta_u = ',F6.4,4X,'delta_v = ',F6.4,4X,'delta_w = ',F6.4)
1609126 FORMAT ('     Upstream-Scheme is used if Upstream-differences fall short', &
1610                  ' of'/                                                       &
1611            '     delta_e = ',F6.4,4X,'delta_pt = ',F6.4)
1612127 FORMAT ('     The following absolute overshoot differences are tolerated:'/&
1613            '     delta_u = ',F6.4,4X,'delta_v = ',F6.4,4X,'delta_w = ',F6.4)
1614128 FORMAT ('     The following absolute overshoot differences are tolerated:'/&
1615            '     delta_e = ',F6.4,4X,'delta_pt = ',F6.4)
1616129 FORMAT (' --> Additional prognostic equation for the specific humidity')
1617130 FORMAT (' --> Additional prognostic equation for the total water content')
1618131 FORMAT (' --> Parameterization of condensation processes via (0%-or100%)')
1619132 FORMAT (' --> Parameterization of long-wave radiation processes via'/ &
1620            '     effective emissivity scheme')
1621133 FORMAT (' --> Precipitation parameterization via Kessler-Scheme')
1622134 FORMAT (' --> Additional prognostic equation for a passive scalar')
1623135 FORMAT (' --> Solve perturbation pressure via multigrid method (', &
1624                  A,'-cycle)'/ &
1625            '     number of grid levels:                   ',I2/ &
1626            '     Gauss-Seidel red/black iterations:       ',I2)
1627136 FORMAT ('     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1628                  I3,')')
1629137 FORMAT ('     level data gathered on PE0 at level:     ',I2/ &
1630            '     gridpoints of coarsest subdomain (x,y,z): (',I3,',',I3,',', &
1631                  I3,')'/ &
1632            '     gridpoints of coarsest domain (x,y,z):    (',I3,',',I3,',', &
1633                  I3,')')
1634138 FORMAT ('     Using hybrid version for 1d-domain-decomposition')
1635139 FORMAT (' --> Loop optimization method: ',A)
1636140 FORMAT ('     maximum residual allowed:                ',E10.3)
1637141 FORMAT ('     fixed number of multigrid cycles:        ',I4)
1638142 FORMAT ('     perturbation pressure is calculated at every Runge-Kutta ', &
1639                  'step')
1640143 FORMAT ('     Euler/upstream scheme is used for the SGS turbulent ', &
1641                  'kinetic energy')
1642150 FORMAT (' --> Volume flow at the right and north boundary will be ', &
1643                  'conserved'/ &
1644            '     using the ',A,' mode')
1645151 FORMAT ('     with u_bulk = ',F7.3,' m/s and v_bulk = ',F7.3,' m/s')
1646152 FORMAT (' --> External pressure gradient directly prescribed by the user:',&
1647           /'     ',2(1X,E12.5),'Pa/m in x/y direction', &
1648           /'     starting from dp_level_b =', F8.3, 'm', A /)
1649153 FORMAT (' --> Large-scale vertical motion is used in the ', &
1650                  'prognostic equation for')
1651154 FORMAT ('     the potential temperature')
1652200 FORMAT (//' Run time and time step information:'/ &
1653             ' ----------------------------------'/)
1654201 FORMAT ( ' Timestep:          variable     maximum value: ',F6.3,' s', &
1655             '    CFL-factor: ',F4.2)
1656202 FORMAT ( ' Timestep:       dt = ',F6.3,' s'/)
1657203 FORMAT ( ' Start time:       ',F9.3,' s'/ &
1658             ' End time:         ',F9.3,' s')
1659204 FORMAT ( A,F9.3,' s')
1660205 FORMAT ( A,F9.3,' s',5X,'restart every',17X,F9.3,' s')
1661206 FORMAT (/' Time reached:     ',F9.3,' s'/ &
1662             ' CPU-time used:    ',F9.3,' s     per timestep:               ', &
1663               '  ',F9.3,' s'/                                                 &
1664             '                                   per second of simulated tim', &
1665               'e: ',F9.3,' s')
1666207 FORMAT ( A/' Coupling start time:',F9.3,' s')
1667250 FORMAT (//' Computational grid and domain size:'/ &
1668              ' ----------------------------------'// &
1669              ' Grid length:      dx =    ',F7.3,' m    dy =    ',F7.3, &
1670              ' m    dz =    ',F7.3,' m'/ &
1671              ' Domain size:       x = ',F10.3,' m     y = ',F10.3, &
1672              ' m  z(u) = ',F10.3,' m'/)
1673252 FORMAT (' dz constant up to ',F10.3,' m (k=',I4,'), then stretched by', &
1674              ' factor: ',F5.3/ &
1675            ' maximum dz not to be exceeded is dz_max = ',F10.3,' m'/)
1676254 FORMAT (' Number of gridpoints (x,y,z):  (0:',I4,', 0:',I4,', 0:',I4,')'/ &
1677            ' Subdomain size (x,y,z):        (  ',I4,',   ',I4,',   ',I4,')'/)
1678255 FORMAT (' Subdomains have equal size')
1679256 FORMAT (' Subdomains at the upper edges of the virtual processor grid ', &
1680              'have smaller sizes'/                                          &
1681            ' Size of smallest subdomain:    (  ',I4,',   ',I4,',   ',I4,')')
1682260 FORMAT (/' The model has a slope in x-direction. Inclination angle: ',F6.2,&
1683             ' degrees')
1684270 FORMAT (//' Topography informations:'/ &
1685              ' -----------------------'// &
1686              1X,'Topography: ',A)
1687271 FORMAT (  ' Building size (x/y/z) in m: ',F5.1,' / ',F5.1,' / ',F5.1/ &
1688              ' Horizontal index bounds (l/r/s/n): ',I4,' / ',I4,' / ',I4, &
1689                ' / ',I4)
1690272 FORMAT (  ' Single quasi-2D street canyon of infinite length in ',A, &
1691              ' direction' / &
1692              ' Canyon height: ', F6.2, 'm, ch = ', I4, '.'      / &
1693              ' Canyon position (',A,'-walls): cxl = ', I4,', cxr = ', I4, '.')
1694278 FORMAT (' Topography grid definition convention:'/ &
1695            ' cell edge (staggered grid points'/  &
1696            ' (u in x-direction, v in y-direction))' /)
1697279 FORMAT (' Topography grid definition convention:'/ &
1698            ' cell center (scalar grid points)' /)
1699280 FORMAT (//' Vegetation canopy (drag) model:'/ &
1700              ' ------------------------------'// &
1701              ' Canopy mode: ', A / &
1702              ' Canopy top: ',I4 / &
1703              ' Leaf drag coefficient: ',F6.2 /)
1704281 FORMAT (/ ' Scalar_exchange_coefficient: ',F6.2 / &
1705              ' Scalar concentration at leaf surfaces in kg/m**3: ',F6.2 /)
1706282 FORMAT (' Predefined constant heatflux at the top of the vegetation: ',F6.2,' K m/s')
1707283 FORMAT (/ ' Characteristic levels of the leaf area density:'// &
1708              ' Height:              ',A,'  m'/ &
1709              ' Leaf area density:   ',A,'  m**2/m**3'/ &
1710              ' Gradient:            ',A,'  m**2/m**4'/ &
1711              ' Gridpoint:           ',A)
1712               
1713300 FORMAT (//' Boundary conditions:'/ &
1714             ' -------------------'// &
1715             '                     p                    uv             ', &
1716             '                   pt'// &
1717             ' B. bound.: ',A/ &
1718             ' T. bound.: ',A)
1719301 FORMAT (/'                     ',A// &
1720             ' B. bound.: ',A/ &
1721             ' T. bound.: ',A)
1722303 FORMAT (/' Bottom surface fluxes are used in diffusion terms at k=1')
1723304 FORMAT (/' Top surface fluxes are used in diffusion terms at k=nzt')
1724305 FORMAT (//'    Prandtl-Layer between bottom surface and first ', &
1725               'computational u,v-level:'// &
1726             '       zp = ',F6.2,' m   z0 = ',F6.4,' m   kappa = ',F4.2/ &
1727             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1728306 FORMAT ('       Predefined constant heatflux:   ',F9.6,' K m/s')
1729307 FORMAT ('       Heatflux has a random normal distribution')
1730308 FORMAT ('       Predefined surface temperature')
1731309 FORMAT ('       Predefined constant salinityflux:   ',F9.6,' psu m/s')
1732310 FORMAT (//'    1D-Model:'// &
1733             '       Rif value range:   ',F6.2,' <= rif <=',F6.2)
1734311 FORMAT ('       Predefined constant humidity flux: ',E10.3,' m/s')
1735312 FORMAT ('       Predefined surface humidity')
1736313 FORMAT ('       Predefined constant scalar flux: ',E10.3,' kg/(m**2 s)')
1737314 FORMAT ('       Predefined scalar value at the surface')
1738315 FORMAT ('       Humidity / scalar flux at top surface is 0.0')
1739316 FORMAT ('       Sensible heatflux and momentum flux from coupled ', &
1740                    'atmosphere model')
1741317 FORMAT (//' Lateral boundaries:'/ &
1742            '       left/right:  ',A/    &
1743            '       north/south: ',A)
1744318 FORMAT (/'       outflow damping layer width: ',I3,' gridpoints with km_', &
1745                    'max =',F5.1,' m**2/s')
1746319 FORMAT ('       turbulence recycling at inflow switched on'/ &
1747            '       width of recycling domain: ',F7.1,' m   grid index: ',I4/ &
1748            '       inflow damping height: ',F6.1,' m   width: ',F6.1,' m')
1749320 FORMAT ('       Predefined constant momentumflux:  u: ',F9.6,' m**2/s**2'/ &
1750            '                                          v: ',F9.6,' m**2/s**2')
1751325 FORMAT (//' List output:'/ &
1752             ' -----------'//  &
1753            '    1D-Profiles:'/    &
1754            '       Output every             ',F8.2,' s')
1755326 FORMAT ('       Time averaged over       ',F8.2,' s'/ &
1756            '       Averaging input every    ',F8.2,' s')
1757330 FORMAT (//' Data output:'/ &
1758             ' -----------'/)
1759331 FORMAT (/'    1D-Profiles:')
1760332 FORMAT (/'       ',A)
1761333 FORMAT ('       Output every             ',F8.2,' s',/ &
1762            '       Time averaged over       ',F8.2,' s'/ &
1763            '       Averaging input every    ',F8.2,' s')
1764334 FORMAT (/'    2D-Arrays',A,':')
1765335 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1766            '       Output every             ',F8.2,' s  ',A/ &
1767            '       Cross sections at ',A1,' = ',A/ &
1768            '       scalar-coordinates:   ',A,' m'/)
1769336 FORMAT (/'    3D-Arrays',A,':')
1770337 FORMAT (/'       Arrays: ',A/ &
1771            '       Output every             ',F8.2,' s  ',A/ &
1772            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1773338 FORMAT ('       Compressed data output'/ &
1774            '       Decimal precision: ',A/)
1775339 FORMAT ('       No output during initial ',F8.2,' s')
1776340 FORMAT (/'    Time series:')
1777341 FORMAT ('       Output every             ',F8.2,' s'/)
1778342 FORMAT (/'       ',A2,'-cross-section  Arrays: ',A/ &
1779            '       Output every             ',F8.2,' s  ',A/ &
1780            '       Time averaged over       ',F8.2,' s'/ &
1781            '       Averaging input every    ',F8.2,' s'/ &
1782            '       Cross sections at ',A1,' = ',A/ &
1783            '       scalar-coordinates:   ',A,' m'/)
1784343 FORMAT (/'       Arrays: ',A/ &
1785            '       Output every             ',F8.2,' s  ',A/ &
1786            '       Time averaged over       ',F8.2,' s'/ &
1787            '       Averaging input every    ',F8.2,' s'/ &
1788            '       Upper output limit at    ',F8.2,' m  (GP ',I4,')'/)
1789344 FORMAT ('       Output format: ',A/)
1790345 FORMAT (/'    Scaling lengths for output locations of all subsequent mask IDs:',/ &
1791            '       mask_scale_x (in x-direction): ',F9.3, ' m',/ &
1792            '       mask_scale_y (in y-direction): ',F9.3, ' m',/ &
1793            '       mask_scale_z (in z-direction): ',F9.3, ' m' )
1794346 FORMAT (/'    Masked data output',A,' for mask ID ',I2, ':')
1795347 FORMAT ('       Variables: ',A/ &
1796            '       Output every             ',F8.2,' s')
1797348 FORMAT ('       Variables: ',A/ &
1798            '       Output every             ',F8.2,' s'/ &
1799            '       Time averaged over       ',F8.2,' s'/ &
1800            '       Averaging input every    ',F8.2,' s')
1801349 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
1802            'mask_scale_',A,' predefined by array mask_',I2.2,'_',A,':'/ &
1803            13('       ',8(F8.2,',')/) )
1804350 FORMAT (/'       Output locations in ',A,'-direction: ', &
1805            'all gridpoints along ',A,'-direction (default).' )
1806351 FORMAT (/'       Output locations in ',A,'-direction in multiples of ', &
1807            'mask_scale_',A,' constructed from array mask_',I2.2,'_',A,'_loop:'/ &
1808            '          loop begin:',F8.2,', end:',F8.2,', stride:',F8.2 )
1809#if defined( __dvrp_graphics )
1810360 FORMAT ('    Plot-Sequence with dvrp-software:'/ &
1811            '       Output every      ',F7.1,' s'/ &
1812            '       Output mode:      ',A/ &
1813            '       Host / User:      ',A,' / ',A/ &
1814            '       Directory:        ',A// &
1815            '       The sequence contains:')
1816361 FORMAT (/'       Isosurface of "',A,'"    Threshold value: ', E12.3/ &
1817            '          Isosurface color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
1818362 FORMAT (/'       Slicer plane ',A/ &
1819            '       Slicer limits: [',F6.2,',',F6.2,']')
1820363 FORMAT (/'       Particles'/ &
1821            '          particle size:  ',F7.2,' m')
1822364 FORMAT ('          particle ',A,' controlled by "',A,'" with interval [', &
1823                       F6.2,',',F6.2,']')
1824365 FORMAT (/'       Groundplate color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)'/ &
1825            '       Superelevation along (x,y,z): (',F4.1,',',F4.1,',',F4.1, &
1826                     ')'/ &
1827            '       Clipping limits: from x = ',F9.1,' m to x = ',F9.1,' m'/ &
1828            '                        from y = ',F9.1,' m to y = ',F9.1,' m')
1829366 FORMAT (/'       Topography color: (',F4.2,',',F4.2,',',F4.2,') (R,G,B)')
1830367 FORMAT ('       Polygon reduction for topography: cluster_size = ', I1)
1831#endif
1832#if defined( __spectra )
1833370 FORMAT ('    Spectra:')
1834371 FORMAT ('       Output every ',F7.1,' s'/)
1835372 FORMAT ('       Arrays:     ', 10(A5,',')/                         &
1836            '       Directions: ', 10(A5,',')/                         &
1837            '       height levels  k = ', 20(I3,',')/                  &
1838            '                          ', 20(I3,',')/                  &
1839            '                          ', 20(I3,',')/                  &
1840            '                          ', 20(I3,',')/                  &
1841            '                          ', 19(I3,','),I3,'.'/           &
1842            '       height levels selected for standard plot:'/        &
1843            '                      k = ', 20(I3,',')/                  &
1844            '                          ', 20(I3,',')/                  &
1845            '                          ', 20(I3,',')/                  &
1846            '                          ', 20(I3,',')/                  &
1847            '                          ', 19(I3,','),I3,'.'/           &
1848            '       Time averaged over ', F7.1, ' s,' /                &
1849            '       Profiles for the time averaging are taken every ', &
1850                    F6.1,' s')
1851#endif
1852400 FORMAT (//' Physical quantities:'/ &
1853              ' -------------------'/)
1854410 FORMAT ('    Angular velocity    :   omega = ',E9.3,' rad/s'/  &
1855            '    Geograph. latitude  :   phi   = ',F4.1,' degr'/   &
1856            '    Coriolis parameter  :   f     = ',F9.6,' 1/s'/    &
1857            '                            f*    = ',F9.6,' 1/s')
1858411 FORMAT (/'    Gravity             :   g     = ',F4.1,' m/s**2')
1859412 FORMAT (/'    Reference density in buoyancy terms: ',F8.3,' kg/m**3')
1860413 FORMAT (/'    Reference temperature in buoyancy terms: ',F8.4,' K')
1861415 FORMAT (/'    Cloud physics parameters:'/ &
1862             '    ------------------------'/)
1863416 FORMAT ('        Surface pressure   :   p_0   = ',F7.2,' hPa'/      &
1864            '        Gas constant       :   R     = ',F5.1,' J/(kg K)'/ &
1865            '        Density of air     :   rho_0 = ',F5.3,' kg/m**3'/  &
1866            '        Specific heat cap. :   c_p   = ',F6.1,' J/(kg K)'/ &
1867            '        Vapourization heat :   L_v   = ',E8.2,' J/kg')
1868420 FORMAT (/'    Characteristic levels of the initial temperature profile:'// &
1869            '       Height:        ',A,'  m'/ &
1870            '       Temperature:   ',A,'  K'/ &
1871            '       Gradient:      ',A,'  K/100m'/ &
1872            '       Gridpoint:     ',A)
1873421 FORMAT (/'    Characteristic levels of the initial humidity profile:'// &
1874            '       Height:      ',A,'  m'/ &
1875            '       Humidity:    ',A,'  kg/kg'/ &
1876            '       Gradient:    ',A,'  (kg/kg)/100m'/ &
1877            '       Gridpoint:   ',A)
1878422 FORMAT (/'    Characteristic levels of the initial scalar profile:'// &
1879            '       Height:                  ',A,'  m'/ &
1880            '       Scalar concentration:    ',A,'  kg/m**3'/ &
1881            '       Gradient:                ',A,'  (kg/m**3)/100m'/ &
1882            '       Gridpoint:               ',A)
1883423 FORMAT (/'    Characteristic levels of the geo. wind component ug:'// &
1884            '       Height:      ',A,'  m'/ &
1885            '       ug:          ',A,'  m/s'/ &
1886            '       Gradient:    ',A,'  1/100s'/ &
1887            '       Gridpoint:   ',A)
1888424 FORMAT (/'    Characteristic levels of the geo. wind component vg:'// &
1889            '       Height:      ',A,'  m'/ &
1890            '       vg:          ',A,'  m/s'/ &
1891            '       Gradient:    ',A,'  1/100s'/ &
1892            '       Gridpoint:   ',A)
1893425 FORMAT (/'    Characteristic levels of the initial salinity profile:'// &
1894            '       Height:     ',A,'  m'/ &
1895            '       Salinity:   ',A,'  psu'/ &
1896            '       Gradient:   ',A,'  psu/100m'/ &
1897            '       Gridpoint:  ',A)
1898426 FORMAT (/'    Characteristic levels of the subsidence/ascent profile:'// &
1899            '       Height:      ',A,'  m'/ &
1900            '       w_subs:      ',A,'  m/s'/ &
1901            '       Gradient:    ',A,'  (m/s)/100m'/ &
1902            '       Gridpoint:   ',A)
1903450 FORMAT (//' LES / Turbulence quantities:'/ &
1904              ' ---------------------------'/)
1905451 FORMAT ('   Diffusion coefficients are constant:'/ &
1906            '   Km = ',F6.2,' m**2/s   Kh = ',F6.2,' m**2/s   Pr = ',F5.2)
1907452 FORMAT ('   Mixing length is limited to the Prandtl mixing lenth.')
1908453 FORMAT ('   Mixing length is limited to ',F4.2,' * z')
1909454 FORMAT ('   TKE is not allowed to fall below ',E9.2,' (m/s)**2')
1910455 FORMAT ('   initial TKE is prescribed as ',E9.2,' (m/s)**2')
1911470 FORMAT (//' Actions during the simulation:'/ &
1912              ' -----------------------------'/)
1913471 FORMAT ('    Disturbance impulse (u,v) every :   ',F6.2,' s'/            &
1914            '    Disturbance amplitude           :     ',F4.2, ' m/s'/       &
1915            '    Lower disturbance level         : ',F8.2,' m (GP ',I4,')'/  &
1916            '    Upper disturbance level         : ',F8.2,' m (GP ',I4,')')
1917472 FORMAT ('    Disturbances continued during the run from i/j =',I4, &
1918                 ' to i/j =',I4)
1919473 FORMAT ('    Disturbances cease as soon as the disturbance energy exceeds',&
1920                 1X,F5.3, ' m**2/s**2')
1921474 FORMAT ('    Random number generator used    : ',A/)
1922475 FORMAT ('    The surface temperature is increased (or decreased, ', &
1923                 'respectively, if'/ &
1924            '    the value is negative) by ',F5.2,' K at the beginning of the',&
1925                 ' 3D-simulation'/)
1926476 FORMAT ('    The surface humidity is increased (or decreased, ',&
1927                 'respectively, if the'/ &
1928            '    value is negative) by ',E8.1,' kg/kg at the beginning of', &
1929                 ' the 3D-simulation'/)
1930477 FORMAT ('    The scalar value is increased at the surface (or decreased, ',&
1931                 'respectively, if the'/ &
1932            '    value is negative) by ',E8.1,' kg/m**3 at the beginning of', &
1933                 ' the 3D-simulation'/)
1934480 FORMAT ('    Particles:'/ &
1935            '    ---------'// &
1936            '       Particle advection is active (switched on at t = ', F7.1, &
1937                    ' s)'/ &
1938            '       Start of new particle generations every  ',F6.1,' s'/ &
1939            '       Boundary conditions: left/right: ', A, ' north/south: ', A/&
1940            '                            bottom:     ', A, ' top:         ', A/&
1941            '       Maximum particle age:                 ',F9.1,' s'/ &
1942            '       Advection stopped at t = ',F9.1,' s'/ &
1943            '       Particles are sorted every ',F9.1,' s'/)
1944481 FORMAT ('       Particles have random start positions'/)
1945482 FORMAT ('          Particles are advected only horizontally'/)
1946483 FORMAT ('       Particles have tails with a maximum of ',I3,' points')
1947484 FORMAT ('            Number of tails of the total domain: ',I10/ &
1948            '            Minimum distance between tailpoints: ',F8.2,' m'/ &
1949            '            Maximum age of the end of the tail:  ',F8.2,' s')
1950485 FORMAT ('       Particle data are written on file every ', F9.1, ' s')
1951486 FORMAT ('       Particle statistics are written on file'/)
1952487 FORMAT ('       Number of particle groups: ',I2/)
1953488 FORMAT ('       SGS velocity components are used for particle advection'/ &
1954            '          minimum timestep for advection: ', F7.5/)
1955489 FORMAT ('       Number of particles simultaneously released at each ', &
1956                    'point: ', I5/)
1957490 FORMAT ('       Particle group ',I2,':'/ &
1958            '          Particle radius: ',E10.3, 'm')
1959491 FORMAT ('          Particle inertia is activated'/ &
1960            '             density_ratio (rho_fluid/rho_particle) = ',F5.3/)
1961492 FORMAT ('          Particles are advected only passively (no inertia)'/)
1962493 FORMAT ('          Boundaries of particle source: x:',F8.1,' - ',F8.1,' m'/&
1963            '                                         y:',F8.1,' - ',F8.1,' m'/&
1964            '                                         z:',F8.1,' - ',F8.1,' m'/&
1965            '          Particle distances:  dx = ',F8.1,' m  dy = ',F8.1, &
1966                       ' m  dz = ',F8.1,' m'/)
1967494 FORMAT ('       Output of particle time series in NetCDF format every ', &
1968                    F8.2,' s'/)
1969495 FORMAT ('       Number of particles in total domain: ',I10/)
1970500 FORMAT (//' 1D-Model parameters:'/                           &
1971              ' -------------------'//                           &
1972            '    Simulation time:                   ',F8.1,' s'/ &
1973            '    Run-controll output every:         ',F8.1,' s'/ &
1974            '    Vertical profile output every:     ',F8.1,' s'/ &
1975            '    Mixing length calculation:         ',A/         &
1976            '    Dissipation calculation:           ',A/)
1977502 FORMAT ('    Damping layer starts from ',F7.1,' m (GP ',I4,')'/)
1978
1979
1980 END SUBROUTINE header
Note: See TracBrowser for help on using the repository browser.