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

Last change on this file since 3246 was 3246, checked in by sward, 6 years ago

Added error handling for wrong input parameters

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