source: palm/trunk/SOURCE/parin.f90 @ 3913

Last change on this file since 3913 was 3885, checked in by kanani, 5 years ago

restructure/add location/debug messages

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