source: palm/tags/release-3.7a/SOURCE/header.f90 @ 3877

Last change on this file since 3877 was 482, checked in by raasch, 14 years ago

New:
---
compare_palm_logs is additionally compiled with mbuild -u (Makefile in trunk/UTIL)

make options (mopts) to be set by configuration file implemented (mrun, mbuild)

humidity=.T. is now usable for runs with topography. wall_humidityflux and
wall_scalarflux are the corresponding new parin arrays.
(check_parameters, init_3d_model, parin)

Large scale vertical motion (subsidence/ascent) can be applied to the
prognostic equation for the potential temperature. (check_parameters, header,
Makefile, modules, parin, prognostic_equations, read_var_list, subsidence,
write_var_list)

A simple method for installing and running palm (with limited features)
has been added. (Makefile, palm_simple_install, palm_simple_run)

Changed:


2d-decomposition is default for Cray-XT machines (init_pegrid)

var_ts is replaced by dots_max (modules,init_3d_model)

Every cloud droplet has now an own weighting factor and can be deleted due to
collisions. Condensation and collision of cloud droplets are adjusted
accordingly. (advec_particles)

Collision efficiency for large cloud droplets has changed according to table of
Rogers and Yau. (collision_efficiency)

Errors:


Bugfix for generating serial jobs (subjob)

Bugfix: index problem concerning gradient_level indices removed (header)

Dimension of array stat in cascade change to prevent type problems with
mpi2 libraries (poisfft_hybrid)

Loop was split to make runs reproducible when using ifort compiler.
(disturb_field)

Bugfix: exchange of ghost points for prho included (time_integration)

Bugfix: calculation of time-averaged surface heatfluxes (sum_up_3d_data)

Bugfix: calculation of precipitation_rate (calc_precipitation)

Bugfix: initial data assignments to some dvrp arrays changed due to error
messages from gfortran compiler (modules)

Bugfix: calculation of cloud droplet velocity (advec_particles)

Bugfix: transfer of particles at south/left edge (advec_particles)

Bugfix: calculation of collision_efficiency (collision_efficiency)

Bugfix: initialisation of var_mod (subsidence)

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