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

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

Give informative message in case initialization_actions has been changed in nested child domain.

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