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

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

New vertical stretching procedure has been introduced

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