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

Last change on this file since 4149 was 4131, checked in by monakurppa, 5 years ago

Several changes in the salsa aerosol module:

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