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

Last change on this file since 3642 was 3637, checked in by knoop, 5 years ago

M Makefile

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