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

Last change on this file since 3182 was 3182, checked in by suehring, 6 years ago

New Inifor features: grid stretching, improved command-interface, support start dates in different formats in both YYYYMMDD and YYYYMMDDHH, Ability to manually control input file prefixes (--radiation-prefix, --soil-preifx, --flow-prefix, --soilmoisture-prefix) for compatiblity with DWD forcast naming scheme; GNU-style short and long option; Prepared output of large-scale forcing profiles (no computation yet); Added preprocessor flag netcdf4 to switch output format between netCDF 3 and 4; Updated netCDF variable names and attributes to comply with PIDS v1.9; Inifor bugfixes: Improved compatibility with older Intel Intel compilers by avoiding implicit array allocation; Added origin_lon/_lat values and correct reference time in dynamic driver global attributes; corresponding PALM changes: adjustments to revised Inifor; variables names in dynamic driver adjusted; enable geostrophic forcing also in offline nested mode; variable names in LES-LES and COSMO offline nesting changed; lateral boundary flags for nesting, in- and outflow conditions renamed

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