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

Last change on this file since 2305 was 2304, checked in by suehring, 7 years ago

Bugfix, enable restarts for child domain

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