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

Last change on this file since 3587 was 3569, checked in by kanani, 5 years ago

Fix for biomet output (ticket:757), merge of uv_exposure into biometeorology_mod

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