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

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

Bugfixes related to spinup and radiation model

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