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

Last change on this file since 2852 was 2849, checked in by Giersch, 6 years ago

Position of d3par namelist in parameter file is unimportant now

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