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

Last change on this file since 4069 was 4022, checked in by suehring, 5 years ago

Synthetic turbulence generator: Revise bias correction of imposed perturbations (correction via volume flow can create instabilities in case the mean volume flow is close to zero); Introduce lower limits in calculation of coefficient matrix, else the calculation may become numerically unstable; Impose perturbations every timestep, even though no new set of perturbations is generated in case dt_stg_call /= dt_3d; Implement a gradual decrease of Reynolds stress and length scales above ABL height (within 1 length scale above ABL depth to 1/10) rather than an abrupt decline; Bugfix in non-nested case: use ABL height for parametrized turbulence; Offline nesting: Rename subroutine for ABL height calculation and make it public; Change top boundary condition for pressure back to Neumann

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