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

Last change on this file since 3695 was 3668, checked in by maronga, 5 years ago

removed most_methods circular and lookup. added improved version of palm_csd

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