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

Last change on this file since 2690 was 2600, checked in by raasch, 6 years ago

small changes concerning r2599, cycle number are now three digits wide

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