source: palm/tags/release-6.0/SOURCE/parin.f90 @ 3636

Last change on this file since 3636 was 3473, checked in by suehring, 5 years ago

Bugfix for previous commit

  • Property svn:keywords set to Id
File size: 48.0 KB
Line 
1!> @file parin.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2018 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: parin.f90 3473 2018-10-30 20:50:15Z raasch $
27! Add virtual measurement module
28!
29! 3472 2018-10-30 20:43:50Z suehring
30! Add indoor model (kanani, srissman, tlang),
31! minor formatting
32!
33! 3467 2018-10-30 19:05:21Z suehring
34! Implementation of a new aerosol module salsa.
35!
36! 3448 2018-10-29 18:14:31Z kanani
37! Add biometeorology
38!
39! 3435 2018-10-26 18:25:44Z gronemeier
40! Add mask_k_over_surface
41!
42! 3421 2018-10-24 18:39:32Z gronemeier
43! Added module for data output at surfaces
44!
45! 3355 2018-10-16 14:03:34Z knoop
46! - offline nesting separated from large-scale forcing module
47! - top boundary condition for pressure in offline nesting changed
48!
49! 3343 2018-10-15 10:38:52Z suehring
50! Introduced reading of date_init to inipar.(Russo)
51! Add extra profiles for chemistry (basit)
52!
53! 3294 2018-10-01 02:37:10Z raasch
54! changes concerning modularization of ocean option
55!
56! 3274 2018-09-24 15:42:55Z knoop
57! Modularization of all bulk cloud physics code components
58!
59! 3248 2018-09-14 09:42:06Z sward
60! Minor formating changes
61!
62! 3246 2018-09-13 15:14:50Z sward
63! Added error handling for input namelist via parin_fail_message
64!
65! 3240 2018-09-12 12:04:40Z Giersch
66! A check that controls the number of user-defined profiles on the restart file
67! with the one given for the current run has been added.
68!
69! 3204 2018-08-23 10:49:05Z raasch
70! additional check for nz
71!
72! 3184 2018-07-27 17:13:09Z suehring
73! Bugfix, comment setting of chemistry and passive scalar top boundary condition
74! in case of offline nesting
75!
76! 3183 2018-07-27 14:25:55Z suehring
77! Rename variables and boundary conditions in mesoscale-offline nesting mode
78!
79! 3182 2018-07-27 13:36:03Z suehring
80! Added multi agent system
81!
82! 3157 2018-07-19 21:08:49Z maronga
83! added use_free_convection_scaling
84!
85! 3083 2018-06-19 14:03:12Z gronemeier
86! Added rans_const_c and rans_const_sigma as input parameters (TG)
87!
88! 3065 2018-06-12 07:03:02Z Giersch
89! New initialization parameters added
90!
91! 3049 2018-05-29 13:52:36Z Giersch
92! Error messages revised
93!
94! 3045 2018-05-28 07:55:41Z Giersch
95! z_max_do2d removed, error messages revised
96!
97! 2995 2018-04-19 12:13:16Z Giersch
98! time_since_reference_point must be calculated/initialized before the first 
99! call of functions related to the radiation model which occur in
100! time_integration_spinup or time_integration
101!
102! 2980 2018-04-17 15:19:27Z suehring
103! Revise message call
104!
105! 2975 2018-04-16 15:22:20Z suehring
106! - Informative message when initializing_actions has been changed
107!   to set_constant_profile in child domain
108! - Change location in message call
109!
110! 2967 2018-04-13 11:22:08Z raasch
111! bugfix: missing parallel cpp-directives added
112!
113! 2941 2018-04-03 11:54:58Z kanani
114! Fix for spinup in case of restart run
115!
116! 2938 2018-03-27 15:52:42Z suehring
117! Change initialization in case child domain should be initialized with Inifor.
118!
119! 2936 2018-03-27 14:49:27Z suehring
120! inipar renamed to initialization_parameters.
121! d3par renamed to runtime_parameters.
122!
123! 2921 2018-03-22 15:05:23Z Giersch
124! Activation of spinup has been moved from lsm/usm_parin to parin itself
125!
126! 2906 2018-03-19 08:56:40Z Giersch
127! ENVIRONMENT variables read/write_svf has been added
128!
129! 2894 2018-03-15 09:17:58Z Giersch
130! read_var_list has been renamed to rrd_global, all module related _parin
131! routines are called before reading the global restart data to overwrite them
132! in case of restart runs
133!
134! 2881 2018-03-13 16:24:40Z suehring
135! Added flag for switching on/off calculation of soil moisture
136!
137! 2849 2018-03-05 10:49:33Z Giersch
138! Position of d3par namelist in parameter file is unimportant now
139!
140! 2826 2018-02-21 12:39:28Z Giersch
141! Bugfix in setting the default boundary conditions for nest domains
142!
143! 2817 2018-02-19 16:32:21Z knoop
144! Preliminary gust module interface implemented
145!
146! 2773 2018-01-30 14:12:54Z suehring
147! Nesting for chemical species implemented
148!
149! 2766 2018-01-22 17:17:47Z kanani
150! Removed preprocessor directive __chem
151!
152! 2718 2018-01-02 08:49:38Z maronga
153! Corrected "Former revisions" section
154!
155! 2696 2017-12-14 17:12:51Z kanani
156! Change in file header (GPL part)
157! Implementation of uv exposure model (FK)
158! Added rans_mode and turbulence_closure to inipar (TG)
159! Implementation of chemistry module
160! Sorting of USE list (FK)
161! Forcing implemented, and initialization with inifor (MS)
162!
163! 2600 2017-11-01 14:11:20Z raasch
164! some comments added and variables renamed concerning r2599
165!
166! 2599 2017-11-01 13:18:45Z hellstea
167! The i/o grouping is updated to work correctly also in nested runs.
168!
169! 2575 2017-10-24 09:57:58Z maronga
170! Renamed phi -> latitude, added longitude
171!
172! 2563 2017-10-19 15:36:10Z Giersch
173! Changed position where restart files are closed.
174!
175! 2550 2017-10-16 17:12:01Z boeske
176! Added complex_terrain
177!
178! 2544 2017-10-13 18:09:32Z maronga
179! Moved day_of_year_init and time_utc_init to inipar.
180!
181! 2397 2017-09-04 16:22:48Z suehring
182! Enable initialization of 3d model by user in the child domain.
183!
184! 2375 2017-08-29 14:10:28Z schwenkel
185! Added aerosol initialization for bulk microphysics
186!
187! 2372 2017-08-25 12:37:32Z sward
188! y_shift added to namelist
189!
190! 2365 2017-08-21 14:59:59Z kanani
191! Vertical grid nesting: add vnest_start_time to d3par (SadiqHuq)
192!
193! 2339 2017-08-07 13:55:26Z gronemeier
194! corrected timestamp in header
195!
196! 2338 2017-08-07 12:15:38Z gronemeier
197! Modularize 1D model
198!
199! 2310 2017-07-11 09:37:02Z gronemeier
200! Bugfix: re-arranged call for error messages for ENVPAR file
201!
202! 2304 2017-07-04 14:35:55Z suehring
203! Bugfix, enable restarts for child domain.
204!
205! 2298 2017-06-29 09:28:18Z raasch
206! -return_addres, return_username in ENVPAR, -cross_ts_uymax, cross_ts_uymin in
207! d3par
208!
209! 2296 2017-06-28 07:53:56Z maronga
210! Added parameters for model spinup
211!
212! 2292 2017-06-20 09:51:42Z schwenkel
213! Implementation of new microphysic scheme: cloud_scheme = 'morrison'
214! includes two more prognostic equations for cloud drop concentration (nc) 
215! and cloud water content (qc).
216!
217! 2267 2017-06-09 09:33:25Z gronemeier
218! Bugfix: removed skipping of reading namelists in case of omitted d3par
219!
220! 2259 2017-06-08 09:09:11Z gronemeier
221! Implemented synthetic turbulence generator
222!
223! 2233 2017-05-30 18:08:54Z suehring
224!
225! 2232 2017-05-30 17:47:52Z suehring
226! typo corrected
227! +wall_salinityflux
228! +tunnel_height, tunnel_lenght, tunnel_width_x, tunnel_width_y,
229!  tunnel_wall_depth
230!
231! 2118 2017-01-17 16:38:49Z raasch
232! -background_communication from inipar
233!
234! 2050 2016-11-08 15:00:55Z gronemeier
235! Implement turbulent outflow condition
236!
237! 2037 2016-10-26 11:15:40Z knoop
238! Anelastic approximation implemented
239!
240! 2035 2016-10-24 15:06:17Z suehring
241! Remove check for npex and npey in nesting case
242!
243! 2011 2016-09-19 17:29:57Z kanani
244! Added flag lsf_exception to allow explicit enabling of large scale forcing
245! together with buildings on flat terrain.
246!
247! 2007 2016-08-24 15:47:17Z kanani
248! Added call to urban surface model for reading of &urban_surface_par
249!
250! 2004 2016-08-24 10:25:59Z suehring
251! Humidity and passive scalar treated separately in nesting mode
252!
253! 2000 2016-08-20 18:09:15Z knoop
254! Forced header and separation lines into 80 columns
255!
256! 1992 2016-08-12 15:14:59Z suehring
257! +top_scalarflux
258!
259! 1960 2016-07-12 16:34:24Z suehring
260! Allocation of s_init
261!
262! 1957 2016-07-07 10:43:48Z suehring
263! flight module added
264!
265! 1955 2016-07-01 12:38:59Z hellstea
266! The parameter intializating_actions is set to 'set_constant_profiles for
267! all nest domains in order to make sure that diagnostic variables are properly
268! initialized for nest domains. Prognostic variables are later initialized by
269! interpolation from the parent domain.
270!
271! 1917 2016-05-27 14:28:12Z witha
272! Initial version of purely vertical nesting introduced.
273!
274! 1914 2016-05-26 14:44:07Z witha
275! Added call to wind turbine model for reading of &wind_turbine_par
276!
277! 1849 2016-04-08 11:33:18Z hoffmann
278! Adapted for modularization of microphysics
279!
280! 1833 2016-04-07 14:23:03Z raasch
281! call of spectra_parin
282!
283! 1831 2016-04-07 13:15:51Z hoffmann
284! turbulence renamed collision_turbulence, drizzle renamed
285! cloud_water_sedimentation
286! curvature_solution_effects removed
287!
288! 1826 2016-04-07 12:01:39Z maronga
289! Added call to radiation model for reading of &radiation_par.
290! Added call to plant canopy model for reading of &canopy_par.
291!
292! 1817 2016-04-06 15:44:20Z maronga
293! Added call to land surface model for reading of &lsm_par
294!
295! 1804 2016-04-05 16:30:18Z maronga
296! Removed code for parameter file check (__check)
297!
298! 1783 2016-03-06 18:36:17Z raasch
299! +netcdf_deflate in d3par, netcdf module and variable names changed
300!
301! 1764 2016-02-28 12:45:19Z raasch
302! cpp-statements for nesting removed, explicit settings of boundary conditions
303! in nest domains,
304! bugfix: npex/npey message moved from inipar to d3par
305! bugfix: check of lateral boundary conditions from check_parameters to here,
306! because they will be already used in init_pegrid and init_grid
307!
308! 1762 2016-02-25 12:31:13Z hellstea
309! Introduction of nested domain feature
310!
311! 1691 2015-10-26 16:17:44Z maronga
312! Added parameter most_method. Renamed prandtl_layer to constant_flux_layer.
313!
314! 1682 2015-10-07 23:56:08Z knoop
315! Code annotations made doxygen readable
316!
317! 1560 2015-03-06 10:48:54Z keck
318! +recycling_yshift
319!
320! 1496 2014-12-02 17:25:50Z maronga
321! Renamed: "radiation -> "cloud_top_radiation"
322!
323! 1484 2014-10-21 10:53:05Z kanani
324! Changes due to new module structure of the plant canopy model:
325!   canopy-model related parameters moved to new package canopy_par in
326!   subroutine package_parin
327!
328! 1429 2014-07-15 12:53:45Z knoop
329! +ensemble_member_nr to prepare the random_generator for ensemble runs
330!
331! 1402 2014-05-09 14:25:13Z raasch
332! location messages modified, progress_bar_disabled included in envpar-NAMELIST
333!
334! 1384 2014-05-02 14:31:06Z raasch
335! location messages added
336!
337! 1365 2014-04-22 15:03:56Z boeske
338! Usage of large scale forcing enabled:
339! +use_subsidence_tendencies
340!
341! 1361 2014-04-16 15:17:48Z hoffmann
342! +call_microphysics_at_all_substeps
343!
344! 1359 2014-04-11 17:15:14Z hoffmann
345! REAL constants provided with KIND-attribute
346!
347! 1353 2014-04-08 15:21:23Z heinze
348! REAL constants provided with KIND-attribute
349!
350! 1327 2014-03-21 11:00:16Z raasch
351! -data_output_format, do3d_compress, do3d_comp_prec
352!
353! 1320 2014-03-20 08:40:49Z raasch
354! ONLY-attribute added to USE-statements,
355! kind-parameters added to all INTEGER and REAL declaration statements,
356! kinds are defined in new module kinds,
357! old module precision_kind is removed,
358! revision history before 2012 removed,
359! comment fields (!:) to be used for variable explanations added to
360! all variable declaration statements
361!
362! 1318 2014-03-17 13:35:16Z raasch
363! +cpu_log_barrierwait in d3par
364!
365! 1301 2014-03-06 13:29:46Z heinze
366! +large_scale_subsidence
367!
368! 1241 2013-10-30 11:36:58Z heinze
369! +nudging
370! +large_scale_forcing
371!
372! 1216 2013-08-26 09:31:42Z raasch
373! +transpose_compute_overlap in inipar
374!
375! 1195 2013-07-01 12:27:57Z heinze
376! Bugfix: allocate ref_state
377!
378! 1179 2013-06-14 05:57:58Z raasch
379! +reference_state in inipar
380!
381! 1159 2013-05-21 11:58:22Z fricke
382! +use_cmax
383!
384! 1128 2013-04-12 06:19:32Z raasch
385! +background_communication in inipar
386!
387! 1115 2013-03-26 18:16:16Z hoffmann
388! unused variables removed
389!
390! 1092 2013-02-02 11:24:22Z raasch
391! unused variables removed
392!
393! 1065 2012-11-22 17:42:36Z hoffmann
394! +nc, c_sedimentation, limiter_sedimentation, turbulence
395! -mu_constant, mu_constant_value
396!
397! 1053 2012-11-13 17:11:03Z hoffmann
398! necessary expansions according to the two new prognostic equations (nr, qr)
399! of the two-moment cloud physics scheme and steering parameters:
400! +*_init, *_surface, *_surface_initial_change, *_vertical_gradient,
401! +*_vertical_gradient_level, surface_waterflux_*,
402! +cloud_scheme, drizzle, mu_constant, mu_constant_value, ventilation_effect
403!
404! 1036 2012-10-22 13:43:42Z raasch
405! code put under GPL (PALM 3.9)
406!
407! 1015 2012-09-27 09:23:24Z raasch
408! -adjust_mixing_length
409!
410! 1003 2012-09-14 14:35:53Z raasch
411! -grid_matching
412!
413! 1001 2012-09-13 14:08:46Z raasch
414! -cut_spline_overshoot, long_filter_factor, overshoot_limit_*, ups_limit_*
415!
416! 996 2012-09-07 10:41:47Z raasch
417! -use_prior_plot1d_parameters
418!
419! 978 2012-08-09 08:28:32Z fricke
420! -km_damp_max, outflow_damping_width
421! +pt_damping_factor, pt_damping_width
422! +z0h_factor
423!
424! 964 2012-07-26 09:14:24Z raasch
425! -cross_normalized_x, cross_normalized_y, cross_xtext, z_max_do1d,
426! z_max_do1d_normalized
427!
428! 940 2012-07-09 14:31:00Z raasch
429! +neutral in inipar
430!
431! 927 2012-06-06 19:15:04Z raasch
432! +masking_method in inipar
433!
434! 824 2012-02-17 09:09:57Z raasch
435! +curvature_solution_effects in inipar
436!
437! 809 2012-01-30 13:32:58Z maronga
438! Bugfix: replaced .AND. and .NOT. with && and ! in the preprocessor directives
439!
440! 807 2012-01-25 11:53:51Z maronga
441! New cpp directive "__check" implemented which is used by check_namelist_files
442!
443! Revision 1.1  1997/07/24 11:22:50  raasch
444! Initial revision
445!
446!
447! Description:
448! ------------
449!> This subroutine reads variables controling the run from the NAMELIST files
450!>
451!> @todo: Revise max_pr_cs (profiles for chemistry)
452!------------------------------------------------------------------------------!
453 SUBROUTINE parin
454 
455
456    USE arrays_3d,                                                             &
457        ONLY:  pt_init, q_init, ref_state, s_init, sa_init,                    &
458               ug, u_init, v_init, vg
459
460    USE biometeorology_mod,                                                    &
461        ONLY:  biom_parin
462
463    USE bulk_cloud_model_mod,                                                  &
464        ONLY:  bcm_parin
465
466    USE chemistry_model_mod,                                                   &
467        ONLY:  chem_parin
468
469    USE chem_modules
470
471    USE control_parameters
472
473    USE cpulog,                                                                &
474        ONLY:  cpu_log_barrierwait
475
476    USE date_and_time_mod,                                                     &
477        ONLY:  date_init, day_of_year_init, time_utc_init
478
479    USE dvrp_variables,                                                        &
480        ONLY:  local_dvrserver_running
481
482    USE flight_mod,                                                            &
483        ONLY:  flight_parin
484
485    USE grid_variables,                                                        &
486        ONLY:  dx, dy
487
488    USE gust_mod,                                                              &
489        ONLY: gust_parin
490
491    USE indices,                                                               &
492        ONLY:  nx, ny, nz
493
494    USE indoor_model_mod,                                                      &
495        ONLY:  im_parin
496
497    USE kinds
498
499    USE land_surface_model_mod,                                                &
500        ONLY: lsm_parin
501
502    USE model_1d_mod,                                                          &
503        ONLY:  damp_level_1d, dt_pr_1d, dt_run_control_1d, end_time_1d
504
505    USE multi_agent_system_mod,                                                &
506        ONLY:  mas_parin
507
508    USE nesting_offl_mod,                                                      &
509        ONLY:  nesting_offl_parin
510       
511    USE netcdf_interface,                                                      &
512        ONLY:  netcdf_data_format, netcdf_deflate, netcdf_precision
513
514    USE ocean_mod,                                                             &
515        ONLY:  ocean_parin
516
517    USE pegrid
518
519    USE plant_canopy_model_mod,                                                &
520         ONLY: pcm_parin
521
522    USE pmc_interface,                                                         &
523        ONLY:  nested_run, nesting_mode
524
525    USE profil_parameter,                                                      &
526        ONLY:  cross_profiles, profile_columns, profile_rows
527
528    USE progress_bar,                                                          &
529        ONLY :  progress_bar_disabled
530
531    USE radiation_model_mod,                                                   &
532        ONLY: radiation_parin
533
534    USE read_restart_data_mod,                                                 &
535        ONLY:  rrd_global 
536             
537    USE salsa_mod,                                                             &
538        ONLY: salsa_parin         
539
540    USE spectra_mod,                                                           &
541        ONLY :  spectra_parin
542
543    USE statistics,                                                            &
544        ONLY:  hom, hom_sum, pr_palm, region, statistic_regions
545
546    USE surface_output_mod,                                                    &
547        ONLY:  surface_output_parin
548
549    USE synthetic_turbulence_generator_mod,                                    &
550        ONLY:  stg_parin
551
552    USE turbulence_closure_mod,                                                &
553        ONLY:  rans_const_c, rans_const_sigma
554
555    USE urban_surface_mod,                                                     &
556        ONLY: usm_parin
557
558    USE uv_exposure_model_mod,                                                 &
559        ONLY:  uvem_parin
560
561    USE vertical_nesting_mod,                                                  &
562        ONLY:  vnest_start_time
563       
564    USE virtual_measurement_mod,                                               &
565        ONLY:  vm_parin
566
567    USE wind_turbine_model_mod,                                                &
568        ONLY:  wtm_parin
569
570
571    IMPLICIT NONE
572
573    CHARACTER (LEN=80) ::  line  !< dummy string that contains the current line of the parameter file
574
575    INTEGER(iwp) ::  global_id      !< process id with respect to MPI_COMM_WORLD
576    INTEGER(iwp) ::  global_procs   !< # of procs with respect to MPI_COMM_WORLD
577    INTEGER(iwp) ::  i              !<
578    INTEGER(iwp) ::  ioerr          !< error flag for open/read/write
579
580    NAMELIST /inipar/  alpha_surface, approximation, bc_e_b,     &
581                       bc_lr, bc_ns, bc_p_b, bc_p_t, bc_pt_b, bc_pt_t, bc_q_b, &
582             bc_q_t,bc_s_b, bc_s_t, bc_uv_b, bc_uv_t,                 &
583             building_height, building_length_x,          &
584             building_length_y, building_wall_left, building_wall_south,       &
585             calc_soil_moisture_during_spinup,                                 &
586             call_psolver_at_all_substeps,  &
587             canyon_height,                                                    &
588             canyon_width_x, canyon_width_y, canyon_wall_left,                 &
589             canyon_wall_south, cfl_factor, cloud_droplets,   &
590             collective_wait, complex_terrain,           &
591             conserve_volume_flow,                                             &
592             conserve_volume_flow_mode, constant_flux_layer,                   &
593             coupling_start_time,             &
594             cycle_mg, damp_level_1d,                                          &
595             data_output_during_spinup,                                        &
596             date_init,                                                        &
597             day_of_year_init,                                                 &
598             dissipation_1d,                                                   &
599             dp_external, dp_level_b, dp_smooth, dpdxy,    &
600             dt, dt_pr_1d, dt_run_control_1d, dt_spinup, dx, dy, dz, dz_max,   &
601             dz_stretch_factor, dz_stretch_level, dz_stretch_level_start,      &
602             dz_stretch_level_end, end_time_1d, ensemble_member_nr, e_init,    &
603             e_min, fft_method, flux_input_mode, flux_output_mode,             &
604             galilei_transformation, humidity,                                 &
605             inflow_damping_height, inflow_damping_width,                      &
606             inflow_disturbance_begin, inflow_disturbance_end,                 &
607             initializing_actions, km_constant,                                &
608             large_scale_forcing, large_scale_subsidence, latitude,            &
609             longitude,                                 &
610             loop_optimization, lsf_exception, masking_method, mg_cycles,      &
611             mg_switch_to_pe0_level, mixing_length_1d, momentum_advec,         &
612             most_method,                                                      &
613             netcdf_precision, neutral, ngsrb,                                 &
614             nsor, nsor_ini, nudging, nx, ny, nz, ocean_mode, omega,           &
615             omega_sor, outflow_source_plane, passive_scalar,                  &
616             prandtl_number, psolver, pt_damping_factor,        &
617             pt_damping_width, pt_reference, pt_surface,                       &
618             pt_surface_initial_change, pt_vertical_gradient,                  &
619             pt_vertical_gradient_level, q_surface, q_surface_initial_change,  &
620             q_vertical_gradient, q_vertical_gradient_level,                   &
621             random_generator, random_heatflux, rans_const_c, rans_const_sigma,&
622             rans_mode,                                                        &
623             rayleigh_damping_factor, rayleigh_damping_height,                 &
624             recycling_width, recycling_yshift,                                &
625             reference_state, residual_limit,                                  &
626             roughness_length,                                                 &
627             scalar_advec,   &
628             scalar_rayleigh_damping,                              &
629             spinup_time, spinup_pt_amplitude, spinup_pt_mean,                 &
630             statistic_regions, subs_vertical_gradient,                        &
631             subs_vertical_gradient_level, surface_heatflux, surface_pressure, &
632             surface_scalarflux, surface_waterflux,                            &
633             s_surface, s_surface_initial_change, s_vertical_gradient,         &
634             s_vertical_gradient_level, time_utc_init, timestep_scheme,        &
635             topography, topography_grid_convention, top_heatflux,             &
636             top_momentumflux_u, top_momentumflux_v,                           &
637             top_scalarflux, transpose_compute_overlap,                        &
638             tunnel_height, tunnel_length, tunnel_width_x, tunnel_width_y,     &
639             tunnel_wall_depth, turbulence_closure,                            &
640             turbulent_inflow, turbulent_outflow,                              &
641             use_subsidence_tendencies, ug_surface, ug_vertical_gradient,      &
642             use_free_convection_scaling,                                      &
643             ug_vertical_gradient_level, use_surface_fluxes, use_cmax,         &
644             use_top_fluxes, use_ug_for_galilei_tr, use_upstream_for_tke,      &
645             uv_heights, u_bulk, u_profile, vg_surface, vg_vertical_gradient,  &
646             vg_vertical_gradient_level, v_bulk, v_profile,&
647             wall_adjustment, wall_heatflux, wall_humidityflux,                &
648             wall_scalarflux, y_shift, zeta_max, zeta_min,  &
649             z0h_factor
650
651    NAMELIST /initialization_parameters/  alpha_surface,                       &
652             approximation, bc_e_b,                                            &
653             bc_lr, bc_ns, bc_p_b, bc_p_t, bc_pt_b, bc_pt_t, bc_q_b,           &
654             bc_q_t,bc_s_b, bc_s_t, bc_uv_b, bc_uv_t,                          &
655             building_height, building_length_x,                               &
656             building_length_y, building_wall_left, building_wall_south,       &
657             calc_soil_moisture_during_spinup,                                 &
658             call_psolver_at_all_substeps,                                     &
659             canyon_height,                                                    &
660             canyon_width_x, canyon_width_y, canyon_wall_left,                 &
661             canyon_wall_south, cfl_factor, cloud_droplets,                    &
662             collective_wait, complex_terrain,                                 &
663             conserve_volume_flow,                                             &
664             conserve_volume_flow_mode, constant_flux_layer,                   &
665             coupling_start_time,                                              &
666             cycle_mg, damp_level_1d,                                          &
667             data_output_during_spinup,                                        &
668             date_init,                                                        &
669             day_of_year_init,                                                 &
670             dissipation_1d,                                                   &
671             dp_external, dp_level_b, dp_smooth, dpdxy,                        &
672             dt, dt_pr_1d, dt_run_control_1d, dt_spinup, dx, dy, dz, dz_max,   &
673             dz_stretch_factor, dz_stretch_level, dz_stretch_level_start,      &
674             dz_stretch_level_end, end_time_1d, ensemble_member_nr, e_init,    &
675             e_min, fft_method, flux_input_mode, flux_output_mode,             &
676             galilei_transformation, humidity,                                 &
677             inflow_damping_height, inflow_damping_width,                      &
678             inflow_disturbance_begin, inflow_disturbance_end,                 &
679             initializing_actions, km_constant,                                &
680             large_scale_forcing, large_scale_subsidence, latitude,            &
681             longitude,                                                        &
682             loop_optimization, lsf_exception, masking_method, mg_cycles,      &
683             mg_switch_to_pe0_level, mixing_length_1d, momentum_advec,         &
684             most_method,                                                      &
685             netcdf_precision, neutral, ngsrb,                                 &
686             nsor, nsor_ini, nudging, nx, ny, nz, ocean_mode, omega,           &
687             omega_sor, outflow_source_plane, passive_scalar,                  &
688             prandtl_number, psolver, pt_damping_factor,                       &
689             pt_damping_width, pt_reference, pt_surface,                       &
690             pt_surface_initial_change, pt_vertical_gradient,                  &
691             pt_vertical_gradient_level, q_surface, q_surface_initial_change,  &
692             q_vertical_gradient, q_vertical_gradient_level,                   &
693             random_generator, random_heatflux, rans_const_c, rans_const_sigma,&
694             rans_mode,                                                        &
695             rayleigh_damping_factor, rayleigh_damping_height,                 &
696             recycling_width, recycling_yshift,                                &
697             reference_state, residual_limit,                                  &
698             roughness_length, scalar_advec,                                   &
699             scalar_rayleigh_damping,                                          &
700             spinup_time, spinup_pt_amplitude, spinup_pt_mean,                 &
701             statistic_regions, subs_vertical_gradient,                        &
702             subs_vertical_gradient_level, surface_heatflux, surface_pressure, &
703             surface_scalarflux, surface_waterflux,                            &
704             s_surface, s_surface_initial_change, s_vertical_gradient,         &
705             s_vertical_gradient_level, time_utc_init, timestep_scheme,        &
706             topography, topography_grid_convention, top_heatflux,             &
707             top_momentumflux_u, top_momentumflux_v,                           &
708             top_scalarflux, transpose_compute_overlap,                        &
709             tunnel_height, tunnel_length, tunnel_width_x, tunnel_width_y,     &
710             tunnel_wall_depth, turbulence_closure,                            &
711             turbulent_inflow, turbulent_outflow,                              &
712             use_subsidence_tendencies, ug_surface, ug_vertical_gradient,      &
713             ug_vertical_gradient_level, use_surface_fluxes, use_cmax,         &
714             use_top_fluxes, use_ug_for_galilei_tr, use_upstream_for_tke,      &
715             use_free_convection_scaling,                                      &
716             uv_heights, u_bulk, u_profile, vg_surface, vg_vertical_gradient,  &
717             vg_vertical_gradient_level, v_bulk, v_profile,                    &
718             wall_adjustment, wall_heatflux, wall_humidityflux,                &
719             wall_scalarflux, y_shift, zeta_max, zeta_min, z0h_factor
720             
721    NAMELIST /d3par/  averaging_interval, averaging_interval_pr,               &
722             cpu_log_barrierwait, create_disturbances,                         &
723             cross_profiles, data_output, data_output_masks,                   &
724             data_output_pr, data_output_2d_on_each_pe, disturbance_amplitude, &
725             disturbance_energy_limit, disturbance_level_b,                    &
726             disturbance_level_t, do2d_at_begin, do3d_at_begin,                &
727             dt, dt_averaging_input, dt_averaging_input_pr,                    &
728             dt_coupling, dt_data_output, dt_data_output_av, dt_disturb,       &
729             dt_domask, dt_dopr, dt_dopr_listing, dt_dots, dt_do2d_xy,         &
730             dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_max, dt_restart,              &
731             dt_run_control,end_time, force_print_header, mask_k_over_surface, &
732             mask_scale_x,                                                     &
733             mask_scale_y, mask_scale_z, mask_x, mask_y, mask_z, mask_x_loop,  &
734             mask_y_loop, mask_z_loop, netcdf_data_format, netcdf_deflate,     &
735             normalizing_region, npex, npey, nz_do3d,                          &
736             profile_columns, profile_rows,     &
737             restart_time, section_xy, section_xz, section_yz,                 &
738             skip_time_data_output, skip_time_data_output_av, skip_time_dopr,  &
739             skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz,          &
740             skip_time_do3d, skip_time_domask, synchronous_exchange,           &
741             termination_time_needed, vnest_start_time
742
743    NAMELIST /runtime_parameters/  averaging_interval, averaging_interval_pr,  &
744             cpu_log_barrierwait, create_disturbances,                         &
745             cross_profiles, data_output, data_output_masks,                   &
746             data_output_pr, data_output_2d_on_each_pe, disturbance_amplitude, &
747             disturbance_energy_limit, disturbance_level_b,                    &
748             disturbance_level_t, do2d_at_begin, do3d_at_begin,                &
749             dt, dt_averaging_input, dt_averaging_input_pr,                    &
750             dt_coupling, dt_data_output, dt_data_output_av, dt_disturb,       &
751             dt_domask, dt_dopr, dt_dopr_listing, dt_dots, dt_do2d_xy,         &
752             dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_max, dt_restart,              &
753             dt_run_control,end_time, force_print_header, mask_k_over_surface, &
754             mask_scale_x,                                                     &
755             mask_scale_y, mask_scale_z, mask_x, mask_y, mask_z, mask_x_loop,  &
756             mask_y_loop, mask_z_loop, netcdf_data_format, netcdf_deflate,     &
757             normalizing_region, npex, npey, nz_do3d,                          &
758             profile_columns, profile_rows,     &
759             restart_time, section_xy, section_xz, section_yz,                 &
760             skip_time_data_output, skip_time_data_output_av, skip_time_dopr,  &
761             skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz,          &
762             skip_time_do3d, skip_time_domask, synchronous_exchange,           &
763             termination_time_needed, vnest_start_time
764
765    NAMELIST /envpar/  progress_bar_disabled, host, local_dvrserver_running,   &
766                       maximum_cpu_time_allowed, maximum_parallel_io_streams,  &
767                       read_svf, revision, run_identifier, tasks_per_node,     &
768                       write_binary, write_svf
769
770!
771!-- First read values of environment variables (this NAMELIST file is
772!-- generated by palmrun)
773    CALL location_message( 'reading environment parameters from ENVPAR', .FALSE. )
774
775    OPEN ( 90, FILE='ENVPAR', STATUS='OLD', FORM='FORMATTED', IOSTAT=ioerr )
776
777    IF ( ioerr /= 0 )  THEN
778       message_string = 'local file ENVPAR not found' //                       &
779                        '&some variables for steering may not be properly set'
780       CALL message( 'parin', 'PA0276', 0, 1, 0, 6, 0 )
781    ELSE
782       READ ( 90, envpar, IOSTAT=ioerr )
783       IF ( ioerr < 0 )  THEN
784          message_string = 'no envpar-NAMELIST found in local file '  //       &
785                           'ENVPAR& or some variables for steering may '  //   &
786                           'not be properly set'
787          CALL message( 'parin', 'PA0278', 0, 1, 0, 6, 0 )
788       ELSEIF ( ioerr > 0 )  THEN
789          message_string = 'errors in local file ENVPAR' //                    &
790                           '&some variables for steering may not be properly set'
791          CALL message( 'parin', 'PA0277', 0, 1, 0, 6, 0 )
792       ENDIF
793       CLOSE ( 90 )
794    ENDIF
795
796    CALL location_message( 'finished', .TRUE. )
797!
798!-- Calculate the number of groups into which parallel I/O is split.
799!-- The default for files which are opened by all PEs (or where each
800!-- PE opens his own independent file) is, that all PEs are doing input/output
801!-- in parallel at the same time. This might cause performance or even more
802!-- severe problems depending on the configuration of the underlying file
803!-- system.
804!-- Calculation of the number of blocks and the I/O group must be based on all
805!-- PEs involved in this run. Since myid and numprocs are related to the
806!-- comm2d communicator, which gives only a subset of all PEs in case of
807!-- nested runs, that information must be inquired again from the global
808!-- communicator.
809!-- First, set the default:
810#if defined( __parallel )
811    CALL MPI_COMM_RANK( MPI_COMM_WORLD, global_id, ierr )
812    CALL MPI_COMM_SIZE( MPI_COMM_WORLD, global_procs, ierr )
813#else
814    global_id    = 0
815    global_procs = 1
816#endif
817    IF ( maximum_parallel_io_streams == -1  .OR.                               &
818         maximum_parallel_io_streams > global_procs )  THEN
819       maximum_parallel_io_streams = global_procs
820    ENDIF
821!
822!-- Now calculate the number of io_blocks and the io_group to which the
823!-- respective PE belongs. I/O of the groups is done in serial, but in parallel
824!-- for all PEs belonging to the same group.
825    io_blocks = global_procs / maximum_parallel_io_streams
826    io_group  = MOD( global_id+1, io_blocks )
827   
828    CALL location_message( 'reading NAMELIST parameters from PARIN', .FALSE. )
829!
830!-- Data is read in parallel by groups of PEs
831    DO  i = 0, io_blocks-1
832       IF ( i == io_group )  THEN
833
834!
835!--       Open the NAMELIST-file which is send with this job
836          CALL check_open( 11 )
837
838!
839!--       Read the control parameters for initialization.
840!--       The namelist "inipar" must be provided in the NAMELIST-file.
841          READ ( 11, initialization_parameters, ERR=10, END=11 )
842          GOTO 14
843         
844 10       BACKSPACE( 11 )
845          READ( 11 , '(A)') line
846          CALL parin_fail_message( 'initialization_parameters', line )
847
848 11       REWIND ( 11 )
849          READ ( 11, inipar, ERR=12, END=13 )
850 
851          message_string = 'namelist inipar is deprecated and will be ' //    &
852                          'removed in near future. & Please use namelist ' // &
853                          'initialization_parameters instead'
854          CALL message( 'parin', 'PA0017', 0, 1, 0, 6, 0 )
855 
856          GOTO 14
857 
858 12       BACKSPACE( 11 )
859          READ( 11 , '(A)') line
860          CALL parin_fail_message( 'inipar', line )
861
862 13       message_string = 'no initialization_parameters-namelist found'
863          CALL message( 'parin', 'PA0272', 1, 2, 0, 6, 0 )
864
865!
866!--       Try to read runtime parameters given by the user for this run
867!--       (namelist "runtime_parameters"). The namelist "runtime_parmeters"   
868!--       can be omitted. In that case default values are used for the         
869!--       parameters.
870 14       line = ' '
871
872          REWIND ( 11 )
873          line = ' '
874          DO WHILE ( INDEX( line, '&runtime_parameters' ) == 0 )
875             READ ( 11, '(A)', END=16 )  line
876          ENDDO
877          BACKSPACE ( 11 )
878
879!
880!--       Read namelist
881          READ ( 11, runtime_parameters, ERR = 15 )
882          GOTO 18
883
884 15       BACKSPACE( 11 )
885          READ( 11 , '(A)') line
886          CALL parin_fail_message( 'runtime_parameters', line )
887
888 16       REWIND ( 11 )
889          line = ' '
890          DO WHILE ( INDEX( line, '&d3par' ) == 0 )
891             READ ( 11, '(A)', END=18 )  line
892          ENDDO
893          BACKSPACE ( 11 )
894
895!
896!--       Read namelist
897          READ ( 11, d3par, ERR = 17, END = 18 )
898
899          message_string = 'namelist d3par is deprecated and will be ' //      &
900                          'removed in near future. &Please use namelist ' //   &
901                          'runtime_parameters instead'
902          CALL message( 'parin', 'PA0487', 0, 1, 0, 6, 0 )
903
904          GOTO 18
905
906 17       BACKSPACE( 11 )
907          READ( 11 , '(A)') line
908          CALL parin_fail_message( 'd3par', line )
909
910 18       CONTINUE
911
912!
913!--       Check for module namelists and read them
914          CALL biom_parin
915          CALL lsm_parin
916          CALL bcm_parin
917          CALL surface_output_parin
918          CALL usm_parin
919          CALL spectra_parin
920          CALL radiation_parin
921          CALL gust_parin
922          CALL mas_parin
923          CALL nesting_offl_parin
924          CALL ocean_parin
925          CALL pcm_parin
926          CALL package_parin
927          CALL wtm_parin
928          CALL flight_parin
929          CALL stg_parin
930          CALL chem_parin
931          CALL uvem_parin
932          CALL im_parin
933          CALL salsa_parin
934          CALL vm_parin
935!
936!--       Read user-defined variables
937          CALL user_parin
938
939!
940!--       If required, read control parameters from restart file (produced by
941!--       a prior run). All PEs are reading from file created by PE0 (see
942!--       check_open)
943          IF ( TRIM( initializing_actions ) == 'read_restart_data' )  THEN
944
945             CALL rrd_global
946!
947!--          Increment the run count
948             runnr = runnr + 1
949!
950!--          In case of a restart run, the number of user-defined profiles on
951!--          the restart file (already stored in max_pr_user) has to match the
952!--          one given for the current run. max_pr_user_tmp is calculated in
953!--          user_parin and max_pr_user is read in via rrd_global.
954             IF ( max_pr_user /= max_pr_user_tmp )  THEN
955                WRITE( message_string, * ) 'the number of user-defined ',      &
956                      'profiles given in data_output_pr (', max_pr_user_tmp,   &
957                      ') does not match the one ',                             &
958                      'found in the restart file (', max_pr_user, ')'
959                CALL message( 'user_parin', 'UI0009', 1, 2, 0, 6, 0 )
960             ENDIF
961          ELSE
962             max_pr_user = max_pr_user_tmp
963          ENDIF
964
965!
966!--       Activate spinup
967          IF ( land_surface .OR. urban_surface )  THEN
968             IF ( spinup_time > 0.0_wp )  THEN
969                coupling_start_time = spinup_time
970                time_since_reference_point = simulated_time - coupling_start_time
971                IF ( spinup_pt_mean == 9999999.9_wp )  THEN
972                   spinup_pt_mean = pt_surface
973                ENDIF
974                end_time = end_time + spinup_time
975                IF ( TRIM( initializing_actions ) /= 'read_restart_data' )     &
976                   spinup = .TRUE.
977             ENDIF
978          ENDIF
979
980!
981!--       In case of nested runs, explicitly set nesting boundary conditions.
982!--       This will overwrite the user settings and basic defaults.
983!--       bc_lr and bc_ns always need to be cyclic for vertical nesting.
984          IF ( nested_run )  THEN
985             IF ( nesting_mode == 'vertical' )  THEN
986                IF (bc_lr /= 'cyclic' .OR. bc_ns /= 'cyclic' )  THEN
987                   WRITE ( message_string, *) 'bc_lr and bc_ns were set to ,', &
988                        'cyclic for vertical nesting'
989                   CALL message( 'parin', 'PA0428', 0, 0, 0, 6, 0 )
990                   bc_lr   = 'cyclic'
991                   bc_ns   = 'cyclic'
992                ENDIF
993                IF ( child_domain )  THEN
994                   bc_uv_t  = 'nested'
995                   bc_pt_t  = 'nested'
996                   bc_q_t   = 'nested'
997                   bc_s_t   = 'nested'
998                   bc_cs_t  = 'nested'
999                   bc_p_t   = 'neumann' 
1000                ENDIF
1001!
1002!--          For other nesting modes only set boundary conditions for
1003!--          nested domains.
1004             ELSE
1005                IF ( child_domain )  THEN
1006                   bc_lr    = 'nested'
1007                   bc_ns    = 'nested'
1008                   bc_uv_t  = 'nested'
1009                   bc_pt_t  = 'nested'
1010                   bc_q_t   = 'nested'
1011                   bc_s_t   = 'nested'
1012                   bc_cs_t  = 'nested'
1013                   bc_p_t   = 'neumann'
1014                ENDIF
1015             ENDIF
1016          ENDIF
1017!
1018!--       Set boundary conditions also in case the model is offline-nested in
1019!--       larger-scale models.
1020          IF ( nesting_offline )  THEN
1021             bc_lr    = 'nesting_offline'
1022             bc_ns    = 'nesting_offline'
1023             bc_uv_t  = 'nesting_offline'
1024             bc_pt_t  = 'nesting_offline'
1025             bc_q_t   = 'nesting_offline'
1026           !  bc_s_t   = 'nesting_offline'  ! scalar boundary condition is not clear
1027           !  bc_cs_t  = 'nesting_offline'  ! same for chemical species
1028!
1029!--          For the pressure set Dirichlet conditions, in contrast to the
1030!--          self nesting. This gives less oscilations within the
1031!--          free atmosphere since the pressure solver has more degrees of
1032!--          freedom. In constrast to the self nesting, this might be
1033!--          justified since the top boundary is far away from the domain
1034!--          of interest.
1035             bc_p_t   = 'dirichlet' !'neumann'
1036          ENDIF
1037
1038!         
1039!--       In case of nested runs, make sure that initializing_actions =
1040!--       'set_constant_profiles' even though the constant-profiles
1041!--       initializations for the prognostic variables will be overwritten
1042!--       by pmci_child_initialize and pmci_parent_initialize. This is,
1043!--       however, important e.g. to make sure that diagnostic variables
1044!--       are set properly. An exception is made in case of restart runs and
1045!--       if user decides to do everything by its own.
1046          IF ( child_domain  .AND.  .NOT. (                                    &
1047               TRIM( initializing_actions ) == 'read_restart_data'      .OR.   &
1048               TRIM( initializing_actions ) == 'set_constant_profiles'  .OR.   &
1049               TRIM( initializing_actions ) == 'by_user' ) )  THEN
1050             message_string = 'initializing_actions = ' //                     &
1051                              TRIM( initializing_actions ) // ' has been ' //  &
1052                              'changed to set_constant_profiles in child ' //  &
1053                              'domain.' 
1054             CALL message( 'parin', 'PA0492', 0, 0, 0, 6, 0 )
1055
1056             initializing_actions = 'set_constant_profiles'
1057          ENDIF           
1058!
1059!--       Check validity of lateral boundary conditions. This has to be done
1060!--       here because they are already used in init_pegrid and init_grid and
1061!--       therefore cannot be check in check_parameters
1062          IF ( bc_lr /= 'cyclic'  .AND.  bc_lr /= 'dirichlet/radiation'  .AND. &
1063               bc_lr /= 'radiation/dirichlet'  .AND.  bc_lr /= 'nested'  .AND. &
1064               bc_lr /= 'nesting_offline' )  THEN
1065             message_string = 'unknown boundary condition: bc_lr = "' // &
1066                              TRIM( bc_lr ) // '"'
1067             CALL message( 'parin', 'PA0049', 1, 2, 0, 6, 0 )
1068          ENDIF
1069          IF ( bc_ns /= 'cyclic'  .AND.  bc_ns /= 'dirichlet/radiation'  .AND. &
1070               bc_ns /= 'radiation/dirichlet'  .AND.  bc_ns /= 'nested'  .AND. &
1071               bc_ns /= 'nesting_offline' )  THEN
1072             message_string = 'unknown boundary condition: bc_ns = "' // &
1073                              TRIM( bc_ns ) // '"'
1074             CALL message( 'parin', 'PA0050', 1, 2, 0, 6, 0 )
1075          ENDIF
1076!
1077!--       Set internal variables used for speed optimization in if clauses
1078          IF ( bc_lr /= 'cyclic' )               bc_lr_cyc    = .FALSE.
1079          IF ( bc_lr == 'dirichlet/radiation' )  bc_lr_dirrad = .TRUE.
1080          IF ( bc_lr == 'radiation/dirichlet' )  bc_lr_raddir = .TRUE.
1081          IF ( bc_ns /= 'cyclic' )               bc_ns_cyc    = .FALSE.
1082          IF ( bc_ns == 'dirichlet/radiation' )  bc_ns_dirrad = .TRUE.
1083          IF ( bc_ns == 'radiation/dirichlet' )  bc_ns_raddir = .TRUE.
1084
1085!
1086!--       Definition of names of areas used for computing statistics. They must
1087!--       be defined at this place, because they are allowed to be redefined by
1088!--       the user in user_parin.
1089          region = 'total domain'
1090
1091!
1092!--       Check in case of initial run, if the grid point numbers are well
1093!--       defined and allocate some arrays which are already needed in
1094!--       init_pegrid or check_parameters. During restart jobs, these arrays
1095!--       will be allocated in rrd_global. All other arrays are allocated
1096!--       in init_3d_model.
1097          IF ( TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
1098
1099             IF ( nx <= 0 )  THEN
1100                WRITE( message_string, * ) 'no value or wrong value given',    &
1101                                           ' for nx: nx=', nx
1102                CALL message( 'parin', 'PA0273', 1, 2, 0, 6, 0 )
1103             ENDIF
1104             IF ( ny <= 0 )  THEN
1105                WRITE( message_string, * ) 'no value or wrong value given',    &
1106                                           ' for ny: ny=', ny
1107                CALL message( 'parin', 'PA0274', 1, 2, 0, 6, 0 )
1108             ENDIF
1109             IF ( nz <= 0 )  THEN
1110                WRITE( message_string, * ) 'no value or wrong value given',    &
1111                                           ' for nz: nz=', nz
1112                CALL message( 'parin', 'PA0275', 1, 2, 0, 6, 0 )
1113             ENDIF
1114
1115!
1116!--          As a side condition, routine flow_statistics require at least 14
1117!--          vertical grid levels (they are used to store time-series data)
1118!>           @todo   Remove this restriction
1119             IF ( nz < 14 )  THEN
1120                WRITE( message_string, * ) 'nz >= 14 is required'
1121                CALL message( 'parin', 'PA0362', 1, 2, 0, 6, 0 )
1122             ENDIF
1123
1124!
1125!--          ATTENTION: in case of changes to the following statement please
1126!--                  also check the allocate statement in routine rrd_global
1127             ALLOCATE( pt_init(0:nz+1), q_init(0:nz+1), s_init(0:nz+1),        &
1128                       ref_state(0:nz+1), sa_init(0:nz+1), ug(0:nz+1),         &
1129                       u_init(0:nz+1), v_init(0:nz+1), vg(0:nz+1),             &
1130                       hom(0:nz+1,2,pr_palm+max_pr_user+max_pr_cs,0:statistic_regions),  &
1131                       hom_sum(0:nz+1,pr_palm+max_pr_user+max_pr_cs,0:statistic_regions) )
1132
1133             hom = 0.0_wp
1134
1135          ENDIF
1136
1137!
1138!--       NAMELIST-file is not needed anymore
1139          CALL close_file( 11 )
1140
1141       ENDIF
1142#if defined( __parallel )
1143       CALL MPI_BARRIER( MPI_COMM_WORLD, ierr )
1144#endif
1145    ENDDO
1146
1147    CALL location_message( 'finished', .TRUE. )
1148
1149 END SUBROUTINE parin
Note: See TracBrowser for help on using the repository browser.