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

Last change on this file since 3421 was 3421, checked in by gronemeier, 5 years ago

new surface-data output; renamed output variables (pt to theta, rho_air to rho, rho_ocean to rho_sea_water)

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