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

Last change on this file since 3542 was 3525, checked in by kanani, 5 years ago

Changes related to clean-up of biometeorology (by dom_dwd_user)

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