source: palm/trunk/SOURCE/modules.f90 @ 861

Last change on this file since 861 was 861, checked in by suehring, 12 years ago

WS5 is available in combination with topography. Version number changed from 3.8 to 3.8a. Modification in init_pt_anomaly.

  • Property svn:keywords set to Id
File size: 67.4 KB
Line 
1 MODULE advection
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6! +wall_flags_0
7! -boundary_flags
8! +nzb_max
9! +adv_sca_1, +adv_mom_1
10!
11! Former revisions:
12! -----------------
13! $Id: modules.f90 861 2012-03-26 14:18:34Z suehring $
14!
15! 849 2012-03-15 10:35:09Z raasch
16! +deleted_particles, deleted_tails, tr.._count_sum, tr.._count_recv_sum in
17! particle_attributes,
18! +de_dx, de_dy, de_dz in arrays_3d,
19! first_call_advec_particles renamed first_call_lpm
20!
21! 828 2012-02-21 12:00:36Z raasch
22! +dissipation_classes, radius_classes, use_kernel_tables,
23! particle feature color renamed class
24!
25! 825 2012-02-19 03:03:44Z raasch
26! +bfactor, curvature_solution_effects, eps_ros, molecular_weight_of_water,
27! vanthoff, -b_cond in cloud_parameters,
28! dopts_num increased to 29, particle attributes speed_x|y|z_sgs renamed
29! rvar1|2|3
30! wang_collision_kernel and turbulence_effects_on_collision replaced by
31! collision_kernel, hall_kernel, palm_kernel, wang_kernel
32!
33! 809 2012-01-30 13:32:58Z marongas
34! Bugfix: replaced .AND. and .NOT. with && and ! in the preprocessor directives
35!
36! 807 2012-01-25 11:53:51Z maronga
37! New cpp directive "__check" implemented which is used by check_namelist_files.
38! New parameter check_restart has been defined which is needed by
39! check_namelist_files only.
40!
41! 805 2012-01-17 15:53:28Z franke
42! Bugfix collective_wait must be out of parallel branch for runs in serial mode
43!
44! 801 2012-01-10 17:30:36Z suehring
45! Dimesion of sums_wsus_ws_l, ! sums_wsvs_ws_l, sums_us2_ws_l, sums_vs2_ws_l,
46! sums_ws2_ws_l, sums_wspts_ws_l, sums_wsqs_ws_l, sums_wssas_ws_l increased.
47! for thread-safe summation in advec_ws.
48!
49! 792 2011-12-01 01:23:23Z raasch
50! particle arrays (particles, parrticles_temp) implemented as pointers,
51! +particles_1, particles_2, sort_count
52!
53! 790 2011-11-29 03:11:20Z raasch
54! +turbulence_effects_on_collision, wang_collision_kernel
55!
56! 785 2011-11-28 09:47:19Z raasch
57! +scalar_rayleigh_damping, rdf_sc
58!
59! 778 2011-11-07 14:18:25Z fricke
60! +gathered_size, subdomain_size
61!
62! 771 2011-10-27 10:56:21Z heinze
63! +lpt_av
64!
65! 767 2011-10-14 06:39:12Z raasch
66! +u_profile, v_profile, uv_heights, use_prescribed_profile_data
67!
68! 759 2011-09-15 13:58:31Z raasch
69! +io_blocks, io_group, maximum_parallel_io_streams,
70! synchronous_exchange moved to control_parameters
71!
72! 743 2011-08-18 16:10:16Z suehring
73! Dimension of sums_wsus_ws_l, sums_wsvs_ws_l, sums_us2_ws_l, sums_vs2_ws_l,
74! sums_ws2_ws_l, sums_wspts_ws_l, sums_wssas_ws_l,sums_wsqs_ws_l needed for
75! statistical evaluation of turbulent fluxes in WS-scheme decreased.
76! 736 2011-08-17 14:13:26Z suehring
77! Dimension of fluxes needed for WS-scheme increased.
78!
79! 722 2011-04-11 06:21:09Z raasch
80! Bugfix: default value for south_border_pe changed to .F.
81!
82! 707 2011-03-29 11:39:40Z raasch
83! +bc_lr_dirrad, bc_lr_raddir, bc_ns_dirrad, bc_ns_raddir, left_border_pe,
84! north_border_pe, right_border_pe, south_border_pe
85! p_sub renamed p_loc
86!
87! 683 2011-02-09 14:25:15Z raasch
88! +synchronous_exchange
89!
90! 673 2011-01-18 16:19:48Z suehring
91! +weight_pres to weight the respective contribution of the Runge-Kutta
92! substeps. +p_sub to buffer the intermediate contributions for Multigrid and
93! SOR.
94!
95! 667 2010-12-23 12:06:00Z suehring/gryschka
96! Removed u_nzb_p1_for_vfc and v_nzb_p1_for_vfc
97! For coupling with different resolution in ocean and atmophere:
98! +nx_a, +nx_o, ny_a, +ny_o, ngp_a, ngp_o, +total_2d_o, +total_2d_a,
99! +coupling_topology
100! Buffer arrays for the left sided advective fluxes added in arrays_3d.
101! +flux_s_u, +flux_s_v, +flux_s_w, +diss_s_u, +diss_s_v, +diss_s_w,
102! +flux_s_pt, +diss_s_pt, +flux_s_e, +diss_s_e, +flux_s_q, +diss_s_q,
103! +flux_s_sa, +diss_s_sa
104! 3d arrays for dissipation control added. (only necessary for vector arch.)
105! +var_x, +var_y, +var_z, +gamma_x, +gamma_y, +gamma_z
106! Default of momentum_advec and scalar_advec changed to 'ws-scheme' .
107! +exchange_mg added in control_parameters to steer the data exchange.
108! Parameters +nbgp, +nxlg, +nxrg, +nysg, +nyng added in indices.
109! flag array +boundary_flags added in indices to steer the degradation of order
110! of the advective fluxes when non-cyclic boundaries are used.
111! MPI-datatypes +type_y, +type_y_int and +type_yz for data_exchange added in
112! pegrid.
113! +sums_wsus_ws_l, +sums_wsvs_ws_l, +sums_us2_ws_l, +sums_vs2_ws_l,
114! +sums_ws2_ws_l, +sums_wspts_ws_l, +sums_wssas_ws_l, +sums_wsqs_ws_l
115! and +weight_substep added in statistics to steer the statistical evaluation
116! of turbulent fluxes in the advection routines.
117! LOGICALS +ws_scheme_sca and +ws_scheme_mom added to get a better performance
118! in prognostic_equations.
119! LOGICAL +dissipation_control control added to steer numerical dissipation
120! in ws-scheme.
121! Changed length of string run_description_header
122!
123! 622 2010-12-10 08:08:13Z raasch
124! +collective_wait in pegrid
125!
126! 600 2010-11-24 16:10:51Z raasch
127! default values of surface_scalarflux and surface_waterflux changed
128! to 9999999.9
129!
130! 580 2010-10-05 13:59:11Z heinze
131! Renaming of ws_vertical_gradient to subs_vertical_gradient,
132! ws_vertical_gradient_level to subs_vertical_gradient_level and
133! ws_vertical_gradient_level_ind to subs_vertical_gradient_level_i
134!
135! 564 2010-09-30 13:18:59Z helmke
136! nc_precision and netcdf_precision dimension changed to 11, all default
137! values of mask_xyz_loop changed to -1.0, dimension of openfile changed to
138! 200+2*max_masks, max_masks changed to 50
139!
140! 553 2010-09-01 14:09:06Z weinreis
141! parameters for masked output are replaced by arrays
142!
143! 531 2010-04-21 06:47:21Z heinze
144! character length of dopr_unit enlarged
145!
146! 519 2010-03-19 05:30:02Z raasch
147! -replace_char, replace_by
148!
149! 493 2010-03-01 08:30:24Z raasch
150! +netcdf_data_format, -netcdf_64bit, -netcdf_64bit_3d, -netcdf_format_mask*,
151! -nc_format_mask, -format_parallel_io
152!
153! 485 2010-02-05 10:57:51Z raasch
154! ngp_3d, ngp_3d_inner changed to 64 bit
155!
156! 449 2010-02-02 11:23:59Z raasch
157! -var_ts: replaced by dots_max,
158! initial data assignments to some dvrp arrays changed due to error messages
159! from gfortran compiler
160! +large_scale_subsidence, ws_vertical_gradient, ws_vertical_gradient_level,
161! ws_vertical_gradient_level_ind, w_subs
162!
163! Branch revisions:
164! -----------------
165! Masked data output: + dt_domask, mask_01~20_x|y|z, mask_01~20_x|y|z_loop,
166! mask_scale|_x|y|z, masks, skip_time_domask
167!
168! control_parameters: scalar and momentum advection set to 'ws-scheme' as
169! default scheme
170! indices: nbgp=3 set as default value
171!
172! Former revisions:
173! -----------------
174! $Id: modules.f90 861 2012-03-26 14:18:34Z suehring $
175!
176! 388 2009-09-23 09:40:33Z raasch
177! +prho, prho_1
178! +bc_lr_cyc, bc_ns_cyc
179! +output_for_t0
180! translation error of actual -> current revisions fixed
181! +q* in dots_label, dots_unit. increased dots_num respectively
182! typographical error in dots_unit fixed
183! +clip_dvrp_*, cluster_size, color_interval, dvrpsize_interval, dvrp_overlap,
184! dvrp_total_overlap, groundplate_color, local_dvrserver_running, n*_dvrp,
185! interval_*_dvrp_prt, isosurface_color, particle_color, particle_dvrpsize,
186! topography color in dvrp_variables,
187! vertical_particle_advection is a 1d-array,
188! +canyon_height, canyon_width_x, canyon_width_y, canyon_wall_left,
189! canyon_wall_south, conserve_volume_flow_mode, coupling_start_time,
190! dp_external, dp_level_b, dp_level_ind_b, dp_smooth, dp_smooth_factor, dpdxy,
191! run_coupled, time_since_reference_point, u_bulk, v_bulk in control_parameters,
192! default value of grid_matching changed to strict
193! +shf_av, qsws_av
194!
195! 217 2008-12-09 18:00:48Z letzel
196! +topography_grid_convention
197! some dvrp-variables changed to single precision, variables for dvrp-mode
198! pathlines added, +target_id, abort_mode, message_string
199!
200! 197 2008-09-16 15:29:03Z raasch
201! allow 100 spectra levels instead of 10 for consistency with
202! define_netcdf_header, +canopy_heat_flux, cthf, lai,
203! +leaf_surface_concentration, scalar_exchange_coefficient, sec, sls
204! +hor_index_bounds, hor_index_bounds_previous_run, id_inflow, id_recycling,
205! inflow_damping_*, mean_inflow_profiles, numprocs_previous_run, nx_on_file,
206! ny_on_file, offset_ocean_*, recycling_plane, recycling_width, turbulent_inflow
207! -myid_char_14
208!
209! 138 2007-11-28 10:03:58Z letzel
210! +drag_coefficient, pch_index, lad_surface, lad_vertical_gradient,
211! lad_vertical_gradient_level, plant_canopy, lad, lad_s, lad_u, lad_v,
212! lad_w, cdc, lad_vertical_gradient_level_ind, canopy_mode
213! +dt_sort_particles, ngp_2dh_s_inner, time_sort_particles, flags,
214! wall_flags_1..10, wall_humidityflux(0:4), wall_qflux(0:4),
215! wall_salinityflux(0:4), wall_scalarflux(0:4)
216!
217! 108 2007-08-24 15:10:38Z letzel
218! +comm_inter, constant_top_momentumflux, coupling_char, coupling_mode,
219! coupling_mode_remote, c_u, c_v, c_w, dt_coupling, e_init, humidity_remote,
220! ngp_xy, nxlu, nysv, port_name, qswst_remote, terminate_coupled,
221! terminate_coupled_remote, time_coupling, top_momentumflux_u|v, type_xy,
222! uswst*, vswst*
223!
224! 97 2007-06-21 08:23:15Z raasch
225! +atmos_ocean_sign, ocean, r, + salinity variables
226! defaults of .._vertical_gradient_levels changed from -1.0 to -9999999.9
227! hydro_press renamed hyp, use_pt_reference renamed use_reference
228!
229! 89 2007-05-25 12:08:31Z raasch
230! +data_output_pr_user, max_pr_user, size of data_output_pr, dopr_index,
231! dopr_initial_index and dopr_unit enlarged,
232! var_hom and var_sum renamed pr_palm
233!
234! 82 2007-04-16 15:40:52Z raasch
235! +return_addres, return_username
236! Cpp-directive lcmuk renamed lc
237!
238! 75 2007-03-22 09:54:05Z raasch
239! +arrays precipitation_amount, precipitation_rate, precipitation_rate_av,
240! rif_wall, z0_av, +arrays u_m_l, u_m_r, etc. for radiation boundary conditions,
241! +loop_optimization, netcdf_64bit_3d, zu_s_inner, zw_w_inner, id_var_zusi_*,
242! id_var_zwwi_*, ts_value, u_nzb_p1_for_vfc, v_nzb_p1_for_vfc, pt_reference,
243! use_pt_reference, precipitation_amount_interval, revision
244! +age_m in particle_type, moisture renamed humidity,
245! -data_output_ts, dots_n, uvmean_outflow, uxrp, vynp,
246! arrays dots_label and dots_unit now dimensioned with dots_max,
247! setting of palm version moved to main program
248!
249! 37 2007-03-01 08:33:54Z raasch
250! +constant_top_heatflux, top_heatflux, use_top_fluxes, +arrays for top fluxes,
251! +nzt_diff, default of bc_pt_t renamed "initial_gradient"
252! Bugfix: p is not a pointer
253!
254! RCS Log replace by Id keyword, revision history cleaned up
255!
256! Revision 1.95  2007/02/11 13:18:30  raasch
257! version 3.1b (last under RCS control)
258!
259! Revision 1.1  1997/07/24 11:21:26  raasch
260! Initial revision
261!
262!
263! Description:
264! ------------
265! Definition of variables for special advection schemes
266!------------------------------------------------------------------------------!
267
268    REAL ::  spl_gamma_x, spl_gamma_y
269
270    REAL, DIMENSION(:), ALLOCATABLE   ::  aex, bex, dex, eex, spl_z_x, spl_z_y
271    REAL, DIMENSION(:,:), ALLOCATABLE ::  spl_tri_x, spl_tri_y, spl_tri_zu, &
272                                          spl_tri_zw
273
274    SAVE
275
276 END MODULE advection
277
278
279
280
281 MODULE array_kind
282
283!------------------------------------------------------------------------------!
284! Description:
285! ------------
286! Definition of type parameters (used for the definition of single or double
287! precision variables)
288!------------------------------------------------------------------------------!
289
290    INTEGER, PARAMETER ::  dpk = SELECTED_REAL_KIND( 12 ), &
291                           spk = SELECTED_REAL_KIND( 6 )
292
293    SAVE
294
295 END MODULE array_kind
296
297
298
299
300 MODULE arrays_3d
301
302!------------------------------------------------------------------------------!
303! Description:
304! ------------
305! Definition of all arrays defined on the computational grid
306!------------------------------------------------------------------------------!
307
308    USE array_kind
309
310    REAL, DIMENSION(:), ALLOCATABLE ::                                         &
311          ddzu, ddzu_pres, dd2zu, dzu, ddzw, dzw, hyp, inflow_damping_factor,  &
312          km_damp_x, km_damp_y, lad, l_grid, pt_init, q_init, rdf, rdf_sc,     &
313          sa_init, ug, u_init, u_nzb_p1_for_vfc, vg, v_init, v_nzb_p1_for_vfc, &
314          w_subs, zu, zw
315
316    REAL, DIMENSION(:,:), ALLOCATABLE ::                                       &
317          c_u, c_v, c_w, dzu_mg, dzw_mg, f1_mg, f2_mg, f3_mg,                  &
318          mean_inflow_profiles, pt_slope_ref, qs, qswst_remote, ts, us, z0,    &
319          flux_s_pt, diss_s_pt, flux_s_e, diss_s_e, flux_s_q, diss_s_q,        &
320          flux_s_sa, diss_s_sa, flux_s_u, flux_s_v, flux_s_w, diss_s_u,        &
321          diss_s_v, diss_s_w, total_2d_o, total_2d_a
322
323    REAL, DIMENSION(:,:), ALLOCATABLE, TARGET ::                               &
324          qsws_1, qsws_2, qswst_1, qswst_2, rif_1, rif_2, saswsb_1, saswst_1,  &
325          shf_1, shf_2, tswst_1, tswst_2, usws_1, usws_2, uswst_1, uswst_2,    &
326          vsws_1, vsws_2, vswst_1, vswst_2
327
328    REAL, DIMENSION(:,:), POINTER ::                                           &
329          qsws, qsws_m, qswst, qswst_m, rif, rif_m, saswsb, saswst, shf,       &
330          shf_m, tswst, tswst_m, usws, uswst, usws_m, uswst_m, vsws, vswst,    &
331          vsws_m, vswst_m
332
333    REAL, DIMENSION(:,:,:), ALLOCATABLE ::                                     &
334          canopy_heat_flux, cdc, d, de_dx, de_dy, de_dz, diss, lad_s, lad_u,   &
335          lad_v, lad_w, lai, l_wall, p_loc, sec, sls, tend, u_m_l, u_m_n,      &
336          u_m_r, u_m_s, v_m_l, v_m_n, v_m_r, v_m_s, w_m_l, w_m_n, w_m_r,       &
337          w_m_s, flux_l_pt, diss_l_pt, flux_l_e, diss_l_e, flux_l_q, diss_l_q, &
338          flux_l_sa, diss_l_sa, flux_l_u, flux_l_v, flux_l_w, diss_l_u,        &
339          diss_l_v, diss_l_w
340
341    REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET ::                             &
342          ql_v, ql_vp
343
344    REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET ::                             &
345          e_1, e_2, e_3, kh_1, kh_2, km_1, km_2, p, prho_1, pt_1, pt_2, pt_3,  &
346          q_1, q_2, q_3, ql_1, ql_2, rho_1, sa_1, sa_2, sa_3, u_1, u_2, u_3,   &
347          v_1, v_2, v_3, vpt_1, vpt_2, w_1, w_2, w_3
348
349    REAL, DIMENSION(:,:,:), POINTER ::                                         &
350          e, e_m, e_p, kh, kh_m, km, km_m, prho, pt, pt_m, pt_p, q, q_m, q_p,  &
351          ql, ql_c, rho, sa, sa_p, te_m, tpt_m, tq_m, tsa_m, tu_m, tv_m, tw_m, &
352          u, u_m, u_p, v, v_m, v_p, vpt, vpt_m, w, w_m, w_p
353
354    REAL, DIMENSION(:,:,:,:), ALLOCATABLE ::  rif_wall
355
356    REAL, DIMENSION(:,:,:), ALLOCATABLE :: var_x, var_y, var_z, gamma_x,       &
357                                           gamma_y, gamma_z
358
359    SAVE
360
361 END MODULE arrays_3d
362
363
364
365
366 MODULE averaging
367
368!------------------------------------------------------------------------------!
369! Description:
370! ------------
371! Definition of variables needed for time-averaging of 2d/3d data
372!------------------------------------------------------------------------------!
373
374    REAL, DIMENSION(:,:), ALLOCATABLE ::  lwp_av, precipitation_rate_av,       &
375                                          qsws_av, shf_av,ts_av, us_av, z0_av
376
377    REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: &
378          e_av, lpt_av, p_av, pc_av, pr_av, pt_av, q_av, ql_av, ql_c_av,     &
379          ql_v_av, ql_vp_av, qv_av, rho_av, s_av, sa_av, u_av, v_av, vpt_av, &
380          w_av
381 END MODULE averaging
382
383
384
385 MODULE cloud_parameters
386
387!------------------------------------------------------------------------------!
388! Description:
389! ------------
390! Definition of variables and constants for cloud physics
391!------------------------------------------------------------------------------!
392
393    LOGICAL ::  curvature_solution_effects = .FALSE.
394
395    REAL  ::  bfactor, cp = 1005.0, diff_coeff_l = 0.23E-4,                    &
396              effective_coll_efficiency,                                       &
397              eps_ros = 1.0E-4,             &  ! accuracy of Rosenbrock method
398              l_d_cp, l_d_r, l_d_rv, l_v = 2.5E+06,                            &
399              mass_of_solute = 1.0E-17,              & ! soluted NaCl in kg
400              molecular_weight_of_solute = 0.05844,  & ! mol. mass NaCl (kg/mol)
401              molecular_weight_of_water = 0.01801528,& ! mol. mass H2O (kg/mol)
402              prec_time_const = 0.001, ql_crit = 0.0005, rho_l = 1.0E3,        &
403              r_d = 287.0, r_v = 461.51, thermal_conductivity_l = 2.43E-2,     &
404              vanthoff = 2.0                        ! van't Hoff factor (NaCl)
405
406    REAL, DIMENSION(:), ALLOCATABLE   ::  pt_d_t, t_d_pt
407
408    REAL, DIMENSION(:,:), ALLOCATABLE ::  precipitation_amount, &
409                                          precipitation_rate
410
411    SAVE
412
413 END MODULE cloud_parameters
414
415
416
417
418 MODULE constants
419
420!------------------------------------------------------------------------------!
421! Description:
422! ------------
423! Definition of general constants
424!------------------------------------------------------------------------------!
425
426    REAL    ::  pi = 3.141592654
427    REAL    ::  adv_mom_1, adv_mom_3, adv_mom_5, adv_sca_1, adv_sca_3, adv_sca_5
428   
429
430    SAVE
431
432 END MODULE constants
433
434
435
436
437 MODULE control_parameters
438
439!------------------------------------------------------------------------------!
440! Description:
441! ------------
442! Definition of parameters for program control
443!------------------------------------------------------------------------------!
444
445    TYPE plot_precision
446       CHARACTER (LEN=6) ::  variable
447       INTEGER           ::  precision
448    END TYPE plot_precision
449
450    TYPE(plot_precision), DIMENSION(100) ::  plot_3d_precision =               &
451                        (/ plot_precision( 'u', 2 ), plot_precision( 'v', 2 ), &
452                           plot_precision( 'w', 2 ), plot_precision( 'p', 5 ), &
453                           plot_precision( 'pt', 2 ),                          &
454                           ( plot_precision( ' ', 1 ), i9 = 1,95 ) /)
455
456    TYPE file_status
457       LOGICAL ::  opened, opened_before
458    END TYPE file_status
459   
460    INTEGER, PARAMETER :: mask_xyz_dimension = 100, max_masks = 50
461
462    TYPE(file_status), DIMENSION(200+2*max_masks) ::                &
463                             openfile = file_status(.FALSE.,.FALSE.)
464
465    CHARACTER (LEN=1)    ::  cycle_mg = 'w', timestep_reason = ' '
466    CHARACTER (LEN=2)    ::  coupling_char = ''
467    CHARACTER (LEN=5)    ::  write_binary = 'false'
468    CHARACTER (LEN=6)    ::  grid_matching = 'strict'
469    CHARACTER (LEN=8)    ::  run_date, run_time
470    CHARACTER (LEN=9)    ::  simulated_time_chr
471    CHARACTER (LEN=11)   ::  topography_grid_convention = ' '
472    CHARACTER (LEN=12)   ::  version = ' ', revision = ' '
473    CHARACTER (LEN=16)   ::  conserve_volume_flow_mode = 'default', &
474                             loop_optimization = 'default', &
475                             momentum_advec = 'ws-scheme', &
476                             psolver = 'poisfft', &
477                             scalar_advec = 'ws-scheme'
478    CHARACTER (LEN=20)   ::  bc_e_b = 'neumann', bc_lr = 'cyclic', &
479                             bc_ns = 'cyclic', bc_p_b = 'neumann', &
480                             bc_p_t = 'dirichlet', bc_pt_b = 'dirichlet', &
481                             bc_pt_t = 'initial_gradient', &
482                             bc_q_b = 'dirichlet', bc_q_t = 'neumann', &
483                             bc_s_b = 'dirichlet', bc_s_t = 'neumann', &
484                             bc_sa_t = 'neumann', &
485                             bc_uv_b = 'dirichlet', bc_uv_t = 'dirichlet', &
486                             canopy_mode = 'block', &
487                             coupling_mode = 'uncoupled', &
488                             coupling_mode_remote = 'uncoupled', &
489                             dissipation_1d = 'as_in_3d_model', &
490                             fft_method = 'system-specific', &
491                             mixing_length_1d = 'as_in_3d_model', &
492                             random_generator = 'numerical-recipes', &
493                             return_addres, return_username, &
494                             timestep_scheme = 'runge-kutta-3'
495    CHARACTER (LEN=40)   ::  avs_data_file, topography = 'flat'
496    CHARACTER (LEN=64)   ::  host = ' '
497    CHARACTER (LEN=80)   ::  log_message, run_identifier
498    CHARACTER (LEN=100)  ::  initializing_actions = ' ' 
499    CHARACTER (LEN=110)  ::  run_description_header
500    CHARACTER (LEN=1000) ::  message_string = ' '
501
502    CHARACTER (LEN=7),  DIMENSION(100) ::  do3d_comp_prec = ' '
503    CHARACTER (LEN=10), DIMENSION(10)  ::  data_output_format = ' '
504    CHARACTER (LEN=11), DIMENSION(100) ::  data_output = ' ',    &
505                                           data_output_user = ' ', doav = ' '
506    CHARACTER (LEN=10), DIMENSION(max_masks,100) ::  &
507         data_output_masks = ' ', data_output_masks_user = ' '
508
509    CHARACTER (LEN=10), DIMENSION(300) ::  data_output_pr = ' '
510    CHARACTER (LEN=10), DIMENSION(200) ::  data_output_pr_user = ' '
511    CHARACTER (LEN=20), DIMENSION(11)  ::  netcdf_precision = ' '
512
513    CHARACTER (LEN=10), DIMENSION(max_masks,0:1,100) ::  domask = ' '
514    CHARACTER (LEN=10), DIMENSION(0:1,100) ::  do2d = ' ', do3d = ' '
515
516    INTEGER ::  abort_mode = 1, average_count_pr = 0, average_count_sp = 0, &
517                average_count_3d = 0, current_timestep_number = 0, &
518                coupling_topology = 0, &
519                dist_range = 0, disturbance_level_ind_b, &
520                disturbance_level_ind_t, doav_n = 0, dopr_n = 0, &
521                dopr_time_count = 0, dopts_time_count = 0, &
522                dosp_time_count = 0, dots_time_count = 0, &
523                do2d_xy_n = 0, do2d_xz_n = 0, do2d_yz_n = 0, do3d_avs_n = 0, &
524                dp_level_ind_b = 0, &
525                dvrp_filecount = 0, dz_stretch_level_index, gamma_mg, &
526                grid_level, ibc_e_b, ibc_p_b, ibc_p_t, ibc_pt_b, ibc_pt_t, &
527                ibc_q_b, ibc_q_t, ibc_sa_t, ibc_uv_b, ibc_uv_t, &
528                inflow_disturbance_begin = -1, inflow_disturbance_end = -1, &
529                intermediate_timestep_count, intermediate_timestep_count_max, &
530                io_group = 0, io_blocks = 1, iran = -1234567, &
531                last_dt_change = 0, masks = 0, maximum_grid_level, &
532                maximum_parallel_io_streams = -1, max_pr_user = 0, &
533                mgcycles = 0, mg_cycles = -1, mg_switch_to_pe0_level = 0, mid, &
534                netcdf_data_format = 2, ngsrb = 2, nsor = 20, &
535                nsor_ini = 100, n_sor, normalizing_region = 0, &
536                nz_do1d, nz_do3d = -9999, outflow_damping_width = -1, &
537                pch_index = 0, prt_time_count = 0, recycling_plane, runnr = 0, &
538                skip_do_avs = 0, terminate_coupled = 0, &
539                terminate_coupled_remote = 0, timestep_count = 0
540
541    INTEGER ::  dist_nxl(0:1), dist_nxr(0:1), dist_nyn(0:1), dist_nys(0:1), &
542                do2d_no(0:1) = 0, do2d_xy_time_count(0:1), &
543                do2d_xz_time_count(0:1), do2d_yz_time_count(0:1), &
544                do3d_no(0:1) = 0, do3d_time_count(0:1), &
545                domask_no(max_masks,0:1) = 0, domask_time_count(max_masks,0:1),&
546                lad_vertical_gradient_level_ind(10) = -9999, &
547                mask_size(max_masks,3) = -1, mask_size_l(max_masks,3) = -1, &
548                mask_start_l(max_masks,3) = -1, &
549                pt_vertical_gradient_level_ind(10) = -9999, &
550                q_vertical_gradient_level_ind(10) = -9999, &
551                sa_vertical_gradient_level_ind(10) = -9999, &
552                section(100,3), section_xy(100) = -9999, &
553                section_xz(100) = -9999, section_yz(100) = -9999, &
554                ug_vertical_gradient_level_ind(10) = -9999, &
555                vg_vertical_gradient_level_ind(10) = -9999, &
556                subs_vertical_gradient_level_i(10) = -9999
557
558#if defined ( __check )
559    INTEGER :: check_restart = 0
560#endif
561
562    INTEGER, DIMENSION(:), ALLOCATABLE ::  grid_level_count
563
564    INTEGER, DIMENSION(:,:), ALLOCATABLE   ::  mask_i, mask_j, mask_k
565    INTEGER, DIMENSION(:,:), ALLOCATABLE   ::  &
566                mask_i_global, mask_j_global, mask_k_global
567
568    LOGICAL ::  adjust_mixing_length = .FALSE., avs_output = .FALSE., &
569                bc_lr_cyc =.TRUE., bc_lr_dirrad = .FALSE., &
570                bc_lr_raddir = .FALSE., bc_ns_cyc = .TRUE., &
571                bc_ns_dirrad = .FALSE., bc_ns_raddir = .FALSE., &
572                call_psolver_at_all_substeps = .TRUE., &
573                cloud_droplets = .FALSE., cloud_physics = .FALSE., &
574                conserve_volume_flow = .FALSE., constant_diffusion = .FALSE., &
575                constant_heatflux = .TRUE., constant_top_heatflux = .TRUE., &
576                constant_top_momentumflux = .FALSE., &
577                constant_top_salinityflux = .TRUE., &
578                constant_waterflux = .TRUE., create_disturbances = .TRUE., &
579                cut_spline_overshoot = .TRUE., &
580                data_output_2d_on_each_pe = .TRUE., & 
581                dissipation_control = .FALSE., disturbance_created = .FALSE., &
582                do2d_at_begin = .FALSE., do3d_at_begin = .FALSE., &
583                do3d_compress = .FALSE., do_sum = .FALSE., &
584                dp_external = .FALSE., dp_smooth = .FALSE., &
585                dt_changed = .FALSE., dt_fixed = .FALSE., &
586                dt_3d_reached, dt_3d_reached_l, exchange_mg = .FALSE., &
587                first_call_lpm = .TRUE., &
588                force_print_header = .FALSE., galilei_transformation = .FALSE.,&
589                humidity = .FALSE., humidity_remote = .FALSE., &
590                inflow_l = .FALSE., inflow_n = .FALSE., &
591                inflow_r = .FALSE., inflow_s = .FALSE., &
592                iso2d_output = .FALSE., large_scale_subsidence = .FALSE., &
593                mg_switch_to_pe0 = .FALSE., &
594                netcdf_output = .FALSE., ocean = .FALSE., &
595                outflow_l = .FALSE., outflow_n = .FALSE., outflow_r = .FALSE., &
596                outflow_s = .FALSE., passive_scalar = .FALSE., &
597                plant_canopy = .FALSE., prandtl_layer = .TRUE., &
598                precipitation = .FALSE., profil_output = .FALSE., &
599                radiation = .FALSE., random_heatflux = .FALSE., &
600                run_control_header = .FALSE., run_coupled = .TRUE., &
601                scalar_rayleigh_damping = .TRUE., sloping_surface = .FALSE., &
602                stop_dt = .FALSE., synchronous_exchange = .FALSE., &
603                terminate_run = .FALSE., turbulent_inflow = .FALSE., &
604                use_prescribed_profile_data = .FALSE., &
605                use_prior_plot1d_parameters = .FALSE., use_reference = .FALSE.,&
606                use_surface_fluxes = .FALSE., use_top_fluxes = .FALSE., &
607                use_ug_for_galilei_tr = .TRUE., use_upstream_for_tke = .FALSE.,&
608                wall_adjustment = .TRUE., ws_scheme_sca = .FALSE., &
609                ws_scheme_mom = .FALSE.
610
611    LOGICAL ::  data_output_xy(0:1) = .FALSE., data_output_xz(0:1) = .FALSE., &
612                data_output_yz(0:1) = .FALSE.
613
614    REAL ::  advected_distance_x = 0.0, advected_distance_y = 0.0, &
615             alpha_surface = 0.0, asselin_filter_factor = 0.1, &
616             atmos_ocean_sign = 1.0, &
617             averaging_interval = 0.0, averaging_interval_pr = 9999999.9, &
618             averaging_interval_sp = 9999999.9, bc_pt_t_val, bc_q_t_val, &
619             bottom_salinityflux = 0.0, &
620             building_height = 50.0, building_length_x = 50.0, &
621             building_length_y = 50.0, building_wall_left = 9999999.9, &
622             building_wall_south = 9999999.9, canyon_height = 50.0, &
623             canyon_width_x = 9999999.9, canyon_width_y = 9999999.9, &
624             canyon_wall_left = 9999999.9, canyon_wall_south = 9999999.9, &
625             cthf = 0.0, cfl_factor = -1.0, cos_alpha_surface, &
626             coupling_start_time, disturbance_amplitude = 0.25, &
627             disturbance_energy_limit = 0.01, &
628             disturbance_level_b = -9999999.9, &
629             disturbance_level_t = -9999999.9, &
630             dp_level_b = 0.0, drag_coefficient = 0.0, &
631             dt = -1.0, dt_averaging_input = 0.0, &
632             dt_averaging_input_pr = 9999999.9, dt_coupling = 9999999.9, &
633             dt_data_output = 9999999.9, &
634             dt_data_output_av = 9999999.9, dt_disturb = 9999999.9, &
635             dt_dopr = 9999999.9, dt_dopr_listing = 9999999.9, &
636             dt_dopts = 9999999.9, dt_dosp = 9999999.9, dt_dots = 9999999.9, &
637             dt_do2d_xy = 9999999.9, dt_do2d_xz = 9999999.9, &
638             dt_do2d_yz = 9999999.9, dt_do3d = 9999999.9, dt_dvrp = 9999999.9, &
639             dt_max = 20.0, dt_restart = 9999999.9, &
640             dt_run_control = 60.0, dt_3d = -1.0, dz = -1.0, &
641             dz_max = 9999999.9, dz_stretch_factor = 1.08, &
642             dz_stretch_level = 100000.0, e_init = 0.0, e_min = 0.0, &
643             end_time = 0.0, &
644             f = 0.0, fs = 0.0, g = 9.81, inflow_damping_height = 9999999.9, &
645             inflow_damping_width = 9999999.9, kappa = 0.4, km_constant = -1.0,&
646             km_damp_max = -1.0, lad_surface = 0.0,  &
647             leaf_surface_concentration = 0.0, long_filter_factor = 0.0, &
648             mask_scale_x = 1.0, mask_scale_y = 1.0, mask_scale_z = 1.0, &
649             maximum_cpu_time_allowed = 0.0,  &
650             molecular_viscosity = 1.461E-5, &
651             old_dt = 1.0E-10, omega = 7.29212E-5, omega_sor = 1.8, &
652             overshoot_limit_e = 0.0, overshoot_limit_pt = 0.0, &
653             overshoot_limit_u = 0.0, overshoot_limit_v = 0.0, &
654             overshoot_limit_w = 0.0, particle_maximum_age = 9999999.9, &
655             phi = 55.0, prandtl_number = 1.0, &
656             precipitation_amount_interval = 9999999.9, prho_reference, &
657             pt_reference = 9999999.9, pt_slope_offset = 0.0, &
658             pt_surface = 300.0, pt_surface_initial_change = 0.0, &
659             q_surface = 0.0, q_surface_initial_change = 0.0, &
660             rayleigh_damping_factor = -1.0, rayleigh_damping_height = -1.0, &
661             recycling_width = 9999999.9, residual_limit = 1.0E-4, &
662             restart_time = 9999999.9, rho_reference, rho_surface, &
663             rif_max = 1.0, rif_min = -5.0, roughness_length = 0.1, &
664             sa_surface = 35.0, scalar_exchange_coefficient = 0.0, &
665             simulated_time = 0.0, simulated_time_at_begin, sin_alpha_surface, &
666             skip_time_data_output = 0.0, skip_time_data_output_av = 9999999.9,&
667             skip_time_dopr = 9999999.9, skip_time_dosp = 9999999.9, &
668             skip_time_do2d_xy = 9999999.9, skip_time_do2d_xz = 9999999.9, &
669             skip_time_do2d_yz = 9999999.9, skip_time_do3d = 9999999.9, &
670             surface_heatflux = 9999999.9, surface_pressure = 1013.25, &
671             surface_scalarflux = 9999999.9, surface_waterflux = 9999999.9, &
672             s_surface = 0.0, s_surface_initial_change = 0.0, &
673             termination_time_needed = -1.0, time_coupling = 0.0, &
674             time_disturb = 0.0, time_dopr = 0.0, time_dopr_av = 0.0, &
675             time_dopr_listing = 0.0, time_dopts = 0.0, time_dosp = 0.0, &
676             time_dosp_av = 0.0, time_dots = 0.0, time_do2d_xy = 0.0, &
677             time_do2d_xz = 0.0, time_do2d_yz = 0.0, time_do3d = 0.0, &
678             time_do_av = 0.0, time_do_sla = 0.0, time_dvrp = 0.0, &
679             time_restart = 9999999.9, time_run_control = 0.0,&
680             time_since_reference_point, top_heatflux = 9999999.9, &
681             top_momentumflux_u = 9999999.9, &
682             top_momentumflux_v = 9999999.9, top_salinityflux = 9999999.9, &
683             ug_surface = 0.0, u_bulk = 0.0, u_gtrans = 0.0, &
684             ups_limit_e = 0.0, ups_limit_pt = 0.0, ups_limit_u = 0.0, &
685             ups_limit_v = 0.0, ups_limit_w = 0.0, vg_surface = 0.0, &
686             v_bulk = 0.0, v_gtrans = 0.0, wall_adjustment_factor = 1.8, &
687             z_max_do1d = -1.0, z_max_do1d_normalized = -1.0, z_max_do2d = -1.0
688
689    REAL ::  do2d_xy_last_time(0:1) = -1.0, do2d_xz_last_time(0:1) = -1.0, &
690             do2d_yz_last_time(0:1) = -1.0, dpdxy(1:2) = 0.0, &
691             dt_domask(max_masks) = 9999999.9, lad_vertical_gradient(10) = 0.0,&
692             lad_vertical_gradient_level(10) = -9999999.9, &
693             mask_scale(3), &
694             pt_vertical_gradient(10) = 0.0, &
695             pt_vertical_gradient_level(10) = -9999999.9, &
696             q_vertical_gradient(10) = 0.0, &
697             q_vertical_gradient_level(10) = -1.0, &
698             s_vertical_gradient(10) = 0.0, &
699             s_vertical_gradient_level(10) = -1.0, &
700             sa_vertical_gradient(10) = 0.0, &
701             sa_vertical_gradient_level(10) = -9999999.9, &
702             skip_time_domask(max_masks) = 9999999.9, threshold(20) = 0.0, &
703             time_domask(max_masks) = 0.0, &
704             tsc(10) = (/ 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 /), &
705             u_profile(100) = 9999999.9, uv_heights(100) = 9999999.9, &
706             v_profile(100) = 9999999.9, &
707             ug_vertical_gradient(10) = 0.0, &
708             ug_vertical_gradient_level(10) = -9999999.9, &
709             vg_vertical_gradient(10) = 0.0, &
710             vg_vertical_gradient_level(10) = -9999999.9, &
711             volume_flow(1:2) = 0.0, volume_flow_area(1:2) = 0.0, &
712             volume_flow_initial(1:2) = 0.0, wall_heatflux(0:4) = 0.0, &
713             wall_humidityflux(0:4) = 0.0, wall_qflux(0:4) = 0.0, &
714             wall_salinityflux(0:4) = 0.0, wall_scalarflux(0:4) = 0.0, &
715             subs_vertical_gradient(10) = 0.0, &
716             subs_vertical_gradient_level(10) = -9999999.9
717
718    REAL, DIMENSION(:), ALLOCATABLE ::  dp_smooth_factor
719
720    REAL, DIMENSION(max_masks,mask_xyz_dimension) :: &
721        mask_x = -1.0, mask_y = -1.0, mask_z = -1.0
722    REAL, DIMENSION(max_masks,3) ::                  &
723        mask_x_loop = -1.0, mask_y_loop = -1.0, mask_z_loop = -1.0
724   
725!
726!--    internal mask arrays ("mask,dimension,selection")
727       REAL, DIMENSION(:,:,:), ALLOCATABLE ::  mask, mask_loop
728
729    SAVE
730
731 END MODULE control_parameters
732
733
734
735
736 MODULE cpulog
737
738!------------------------------------------------------------------------------!
739! Description:
740! ------------
741! Definition of variables for cpu-time measurements
742!------------------------------------------------------------------------------!
743
744    REAL ::  initial_wallclock_time
745
746    TYPE logpoint
747       REAL               ::  isum, ivect, mean, mtime, mtimevec, sum, vector
748       INTEGER            ::  counts
749       CHARACTER (LEN=20) ::  place
750    END TYPE logpoint
751
752    TYPE(logpoint), DIMENSION(100) ::  log_point = logpoint( 0.0, 0.0, 0.0,   &
753                                       0.0, 0.0, 0.0, 0.0, 0, ' ' ),          &
754                                       log_point_s = logpoint( 0.0, 0.0, 0.0, &
755                                       0.0, 0.0, 0.0, 0.0, 0, ' ' )
756
757    SAVE
758
759 END MODULE cpulog
760
761
762
763
764 MODULE dvrp_variables
765
766!------------------------------------------------------------------------------!
767! Description:
768! ------------
769! Definition of variables used with dvrp-software
770!------------------------------------------------------------------------------!
771
772    CHARACTER (LEN=10) ::  dvrp_output = 'rtsp', particle_color = 'none', &
773                           particle_dvrpsize = 'none'
774
775    CHARACTER (LEN=20), DIMENSION(10) ::  mode_dvrp = &
776                                     (/ ( '                    ', i9 = 1,10 ) /)
777
778    CHARACTER (LEN=80) ::  dvrp_directory = 'default',                    &
779                           dvrp_file      = 'default',                    &
780                           dvrp_host      = 'origin.rvs.uni-hannover.de', &
781                           dvrp_password  = '********',                   &
782                           dvrp_username  = ' '
783
784    INTEGER ::  cluster_size = 1, dvrp_colortable_entries = 4,                 &
785                dvrp_colortable_entries_prt = 22, islice_dvrp,                 &
786                nx_dvrp, nxl_dvrp, nxr_dvrp, ny_dvrp, nyn_dvrp, nys_dvrp,      &
787                nz_dvrp, pathlines_fadeintime = 5, pathlines_fadeouttime = 5,  &
788                pathlines_linecount = 1000, pathlines_maxhistory = 40,         &
789                pathlines_wavecount = 10, pathlines_wavetime = 50,             &
790                vc_gradient_normals = 0, vc_mode = 0, vc_size_x = 2,           &
791                vc_size_y = 2, vc_size_z = 2
792
793    INTEGER, DIMENSION(10) ::  slicer_position_dvrp
794
795    LOGICAL ::  cyclic_dvrp = .FALSE., dvrp_overlap, dvrp_total_overlap, &
796                local_dvrserver_running, lock_steering_update = .FALSE., &
797                use_seperate_pe_for_dvrp_output = .FALSE.
798
799    REAL    ::  clip_dvrp_l = 9999999.9, clip_dvrp_n = 9999999.9, &
800                clip_dvrp_r = 9999999.9, clip_dvrp_s = 9999999.9, &
801                superelevation = 1.0, superelevation_x = 1.0,     &
802                superelevation_y = 1.0, vc_alpha = 38.0
803
804    REAL, DIMENSION(2) ::  color_interval = (/ 0.0, 1.0 /), &
805                           dvrpsize_interval = (/ 0.0, 1.0 /)
806
807    REAL, DIMENSION(3) ::  groundplate_color = (/ 0.0, 0.6, 0.0 /), &
808                           topography_color = (/ 0.8, 0.7, 0.6 /)
809
810#if defined( __decalpha )
811    REAL, DIMENSION(2,10)  ::  slicer_range_limits_dvrp = RESHAPE( (/       &
812                                -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, &
813                                -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, &
814                                -1.0, 1.0, -1.0, 1.0 /), (/ 2, 10 /) )
815
816    REAL, DIMENSION(3,10)  ::  isosurface_color = RESHAPE( (/                  &
817                                0.9, 0.9, 0.9,  0.8, 0.1, 0.1,  0.1, 0.1, 0.8, &
818                                0.1, 0.8, 0.1,  0.6, 0.1, 0.1,  0.1, 0.1, 0.6, &
819                                0.1, 0.6, 0.1,  0.4, 0.1, 0.1,  0.1, 0.1, 0.4, &
820                                0.1, 0.4, 0.1 /), (/ 3, 10 /) )
821
822    REAL(4), DIMENSION(2,100) ::  interval_values_dvrp, interval_h_dvrp =      &
823                                  RESHAPE( (/ 270.0, 225.0, 225.0, 180.0,      &
824                                               70.0,  25.0,  25.0, -25.0,      &
825                                              ( 0.0, i9 = 1, 192 ) /),         &
826                                           (/ 2, 100 /) ),                     &
827                                  interval_l_dvrp = 0.5, interval_s_dvrp = 1.0,&
828                                  interval_a_dvrp = 0.0,                       &
829                                  interval_values_dvrp_prt,                    &
830                                  interval_h_dvrp_prt = RESHAPE(               &
831                                  (/ 270.0, 225.0, 225.0, 180.0, 70.0, 25.0,   &
832                                     25.0, -25.0, ( 0.0, i9 = 1, 192 ) /),     &
833                                                  (/ 2, 100 /) ),              &
834                                  interval_l_dvrp_prt = 0.5,                   &
835                                  interval_s_dvrp_prt = 1.0,                   &
836                                  interval_a_dvrp_prt = 0.0
837#else
838    REAL, DIMENSION(2,10)     ::  slicer_range_limits_dvrp
839
840    REAL, DIMENSION(3,10)     ::  isosurface_color
841
842    REAL(4), DIMENSION(2,100) ::  interval_values_dvrp,                       &
843                                  interval_values_dvrp_prt, interval_h_dvrp,  &
844                                  interval_h_dvrp_prt, interval_l_dvrp = 0.5, &
845                                  interval_l_dvrp_prt = 0.5, interval_s_dvrp = 1.0, &
846                                  interval_s_dvrp_prt = 1.0, interval_a_dvrp = 0.0, &
847                                  interval_a_dvrp_prt = 0.0
848
849    DATA  slicer_range_limits_dvrp / -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, &
850                                     -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, &
851                                     -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, &
852                                     -1.0, 1.0 /
853
854    DATA  isosurface_color / 0.9, 0.9, 0.9,  0.8, 0.1, 0.1,  0.1, 0.1, 0.8, &
855                             0.1, 0.8, 0.1,  0.6, 0.1, 0.1,  0.1, 0.1, 0.6, &
856                             0.1, 0.6, 0.1,  0.4, 0.1, 0.1,  0.1, 0.1, 0.4, &
857                             0.1, 0.4, 0.1 /
858
859    DATA  interval_h_dvrp / 270.0, 225.0, 225.0, 180.0, 70.0, 25.0, &
860                            25.0, -25.0, 192 * 0.0 /
861
862    DATA  interval_h_dvrp_prt / 270.0, 225.0, 225.0, 180.0, 70.0, 25.0, &
863                                25.0, -25.0, 192 * 0.0 /
864#endif
865
866    REAL(4), DIMENSION(:), ALLOCATABLE ::  xcoor_dvrp, ycoor_dvrp, zcoor_dvrp
867
868    TYPE steering
869       CHARACTER (LEN=20) ::  name
870       REAL(4)            ::  min, max
871       INTEGER            ::  imin, imax
872    END TYPE steering
873
874    TYPE(steering), DIMENSION(:), ALLOCATABLE ::  steering_dvrp
875
876    SAVE
877
878 END MODULE dvrp_variables
879
880
881
882
883 MODULE grid_variables
884
885!------------------------------------------------------------------------------!
886! Description:
887! ------------
888! Definition of grid spacings
889!------------------------------------------------------------------------------!
890
891    REAL ::  ddx, ddx2, dx = 1.0, dx2, ddy, ddy2, dy = 1.0, dy2
892
893    REAL, DIMENSION(:), ALLOCATABLE ::  ddx2_mg, ddy2_mg
894
895    REAL, DIMENSION(:,:), ALLOCATABLE ::  fwxm, fwxp, fwym, fwyp, fxm, fxp,   &
896                                          fym, fyp, wall_e_x, wall_e_y,       &
897                                          wall_u, wall_v, wall_w_x, wall_w_y, &
898                                          zu_s_inner, zw_w_inner
899
900    SAVE
901
902 END MODULE grid_variables
903
904
905
906
907 MODULE indices
908
909!------------------------------------------------------------------------------!
910! Description:
911! ------------
912! Definition of array bounds and number of gridpoints
913!------------------------------------------------------------------------------!
914
915    INTEGER ::  nbgp = 3, ngp_sums, nnx, nx = 0, nx_a, nx_o, nxa, nxl, nxlg,   &
916                nxlu, nxr, nxra, nxrg, nx_on_file, nny, ny = 0, ny_a, ny_o,    &
917                nya, nyn, nyna, nyng, nys, nysg, nysv, ny_on_file, nnz, nz = 0,&
918                nza, nzb, nzb_diff, nzb_max, nzt, nzta, nzt_diff
919
920
921    INTEGER( KIND = SELECTED_INT_KIND(18) ), DIMENSION(:), ALLOCATABLE ::      &
922                ngp_3d, ngp_3d_inner   ! need to have 64 bit for grids > 2E9
923
924    INTEGER, DIMENSION(:), ALLOCATABLE ::                                      &
925                ngp_2dh, nnx_pe, nny_pe, nxl_mg, nxr_mg, nyn_mg, nys_mg, nzt_mg
926
927
928    INTEGER, DIMENSION(:,:), ALLOCATABLE :: ngp_2dh_outer, ngp_2dh_s_inner,    &
929                mg_loc_ind, nzb_diff_s_inner, nzb_diff_s_outer, nzb_diff_u,    &
930                nzb_diff_v, nzb_inner, nzb_outer, nzb_s_inner, nzb_s_outer,    &
931                nzb_u_inner, nzb_u_outer, nzb_v_inner, nzb_v_outer,            &
932                nzb_w_inner, nzb_w_outer, nzb_2d
933
934    INTEGER, DIMENSION(:,:,:), POINTER ::  flags
935
936    INTEGER( KIND = SELECTED_INT_KIND(11) ), DIMENSION(:,:,:), ALLOCATABLE ::  &
937                wall_flags_0  ! need to have 34 Bit
938
939    INTEGER, DIMENSION(:,:,:), ALLOCATABLE,  TARGET ::                         &
940                wall_flags_1, wall_flags_2, wall_flags_3, wall_flags_4,        &
941                wall_flags_5, wall_flags_6, wall_flags_7, wall_flags_8,        &
942                wall_flags_9, wall_flags_10
943
944
945    SAVE
946
947 END MODULE indices
948
949
950
951
952 MODULE interfaces
953
954!------------------------------------------------------------------------------!
955! Description:
956! ------------
957! Interfaces for special subroutines which use optional parameters
958!------------------------------------------------------------------------------!
959
960    INTERFACE
961
962       SUBROUTINE cpu_log( log_event, place, modus, barrierwait )
963
964          USE cpulog
965
966          CHARACTER (LEN=*)           ::  modus, place
967          CHARACTER (LEN=*), OPTIONAL ::  barrierwait
968          TYPE(logpoint)              ::  log_event
969
970       END SUBROUTINE cpu_log
971
972    END INTERFACE
973
974
975
976    INTERFACE
977
978       SUBROUTINE global_min_max ( i1, i2, j1, j2, k1, k2, feld, mode, wert, &
979                                   wert_ijk, wert1, wert1_ijk )
980          CHARACTER (LEN=*), INTENT(IN) ::  mode
981          INTEGER, INTENT(IN)           ::  i1, i2, j1, j2, k1, k2
982          INTEGER                       ::  wert_ijk(3)
983          INTEGER, OPTIONAL             ::  wert1_ijk(3)
984          REAL                          ::  wert
985          REAL, OPTIONAL                ::  wert1
986          REAL, INTENT(IN)              ::  feld(i1:i2,j1:j2,k1:k2)
987
988       END SUBROUTINE global_min_max
989
990    END INTERFACE
991
992    SAVE
993
994 END MODULE interfaces
995
996
997
998 MODULE pointer_interfaces
999
1000!------------------------------------------------------------------------------!
1001! Description:
1002! ------------
1003! Interfaces for subroutines with pointer arguments called in
1004! prognostic_equations
1005!------------------------------------------------------------------------------!
1006
1007    INTERFACE
1008
1009       SUBROUTINE advec_s_bc( sk, sk_char )
1010
1011          CHARACTER (LEN=*), INTENT(IN)   ::  sk_char
1012          REAL, DIMENSION(:,:,:), POINTER ::  sk
1013
1014       END SUBROUTINE advec_s_bc
1015
1016    END INTERFACE
1017
1018
1019    SAVE
1020
1021 END MODULE pointer_interfaces
1022
1023
1024
1025
1026 MODULE model_1d
1027
1028!------------------------------------------------------------------------------!
1029! Description:
1030! ------------
1031! Definition of variables for the 1D-model
1032!------------------------------------------------------------------------------!
1033
1034    INTEGER ::  current_timestep_number_1d = 0, damp_level_ind_1d, &
1035                last_dt_change_1d = 0 
1036
1037    LOGICAL ::  run_control_header_1d = .FALSE., stop_dt_1d = .FALSE. 
1038
1039    REAL ::     damp_level_1d = -1.0, dt_1d = 60.0, dt_max_1d = 300.0, &
1040                dt_pr_1d = 9999999.9, dt_run_control_1d = 60.0, &
1041                end_time_1d = 864000.0, old_dt_1d = 1.0E-10, &
1042                qs1d, simulated_time_1d = 0.0, time_pr_1d = 0.0, &
1043                time_run_control_1d = 0.0, ts1d, us1d, usws1d, usws1d_m, &
1044                vsws1d, vsws1d_m, z01d
1045
1046
1047    REAL, DIMENSION(:), ALLOCATABLE ::  e1d, e1d_m, e1d_p, kh1d, kh1d_m, km1d, &
1048                                        km1d_m, l_black, l1d, l1d_m, rif1d,    &
1049                                        te_e, te_em, te_u, te_um, te_v, te_vm, &
1050                                        u1d, u1d_m, u1d_p, v1d, v1d_m, v1d_p
1051
1052    SAVE
1053
1054 END MODULE model_1d
1055
1056
1057
1058
1059 MODULE netcdf_control
1060
1061!------------------------------------------------------------------------------!
1062! Description:
1063! ------------
1064! Definition of parameters and variables for netcdf control.
1065!------------------------------------------------------------------------------!
1066
1067    USE control_parameters
1068#if defined( __netcdf )
1069    USE netcdf
1070#endif
1071
1072    INTEGER, PARAMETER ::  dopr_norm_num = 7, dopts_num = 29, dots_max = 100
1073
1074    INTEGER ::  dots_num = 23
1075
1076    CHARACTER (LEN=6), DIMENSION(dopr_norm_num) ::  dopr_norm_names =   &
1077         (/ 'wpt0  ', 'ws2   ', 'tsw2  ', 'ws3   ', 'ws2tsw', 'wstsw2', &
1078            'z_i   ' /)
1079
1080    CHARACTER (LEN=6), DIMENSION(dopr_norm_num) ::  dopr_norm_longnames = &
1081         (/ 'wpt0  ', 'w*2   ', 't*w2  ', 'w*3   ', 'w*2t*w', 'w*t*w2',   &
1082            'z_i   ' /)
1083
1084    CHARACTER (LEN=7), DIMENSION(dopts_num) :: dopts_label = &
1085          (/ 'tnpt   ', 'x_     ', 'y_     ', 'z_     ', 'z_abs  ', 'u      ', &
1086             'v      ', 'w      ', 'u"     ', 'v"     ', 'w"     ', 'npt_up ', &
1087             'w_up   ', 'w_down ', 'radius ', 'r_min  ', 'r_max  ', 'npt_max', &
1088             'npt_min', 'x*2    ', 'y*2    ', 'z*2    ', 'u*2    ', 'v*2    ', &
1089             'w*2    ', 'u"2    ', 'v"2    ', 'w"2    ', 'npt*2  ' /)
1090
1091    CHARACTER (LEN=7), DIMENSION(dopts_num) :: dopts_unit = &
1092          (/ 'number ', 'm      ', 'm      ', 'm      ', 'm      ', 'm/s    ', &
1093             'm/s    ', 'm/s    ', 'm/s    ', 'm/s    ', 'm/s    ', 'number ', &
1094             'm/s    ', 'm/s    ', 'm      ', 'm      ', 'm      ', 'number ', &
1095             'number ', 'm2     ', 'm2     ', 'm2     ', 'm2/s2  ', 'm2/s2  ', &
1096             'm2/s2  ', 'm2/s2  ', 'm2/s2  ', 'm2/s2  ', 'number2' /)
1097
1098    CHARACTER (LEN=7), DIMENSION(dots_max) :: dots_label = &
1099          (/ 'E      ', 'E*     ', 'dt     ', 'u*     ', 'th*    ', 'umax   ', &
1100             'vmax   ', 'wmax   ', 'div_new', 'div_old', 'z_i_wpt', 'z_i_pt ', &
1101             'w*     ', 'w"pt"0 ', 'w"pt"  ', 'wpt    ', 'pt(0)  ', 'pt(zp) ', &
1102             'w"u"0  ', 'w"v"0  ', 'w"q"0  ', 'mo_L   ', 'q*     ',            &
1103             ( 'unknown', i9 = 1, dots_max-23 ) /)
1104
1105    CHARACTER (LEN=7), DIMENSION(dots_max) :: dots_unit = &
1106          (/ 'm2/s2  ', 'm2/s2  ', 's      ', 'm/s    ', 'K      ', 'm/s    ', &
1107             'm/s    ', 'm/s    ', 's-1    ', 's-1    ', 'm      ', 'm      ', &
1108             'm/s    ', 'K m/s  ', 'K m/s  ', 'K m/s  ', 'K      ', 'K      ', &
1109             'm2/s2  ', 'm2/s2  ', 'kg m/s ', 'm      ', 'kg/kg  ',            &
1110             ( 'unknown', i9 = 1, dots_max-23 ) /)
1111
1112    CHARACTER (LEN=9), DIMENSION(300) ::  dopr_unit = 'unknown'
1113
1114    CHARACTER (LEN=7), DIMENSION(0:1,100) ::  do2d_unit, do3d_unit
1115
1116    CHARACTER (LEN=16), DIMENSION(25) ::  prt_var_names = &
1117          (/ 'pt_age          ', 'pt_dvrp_size    ', 'pt_origin_x     ', &
1118             'pt_origin_y     ', 'pt_origin_z     ', 'pt_radius       ', &
1119             'pt_speed_x      ', 'pt_speed_y      ', 'pt_speed_z      ', &
1120             'pt_weight_factor', 'pt_x            ', 'pt_y            ', &
1121             'pt_z            ', 'pt_color        ', 'pt_group        ', &
1122             'pt_tailpoints   ', 'pt_tail_id      ', 'pt_density_ratio', &
1123             'pt_exp_arg      ', 'pt_exp_term     ', 'not_used        ', &
1124             'not_used        ', 'not_used        ', 'not_used        ', &
1125             'not_used        ' /)
1126
1127    CHARACTER (LEN=16), DIMENSION(25) ::  prt_var_units = &
1128          (/ 'seconds         ', 'meters          ', 'meters          ', &
1129             'meters          ', 'meters          ', 'meters          ', &
1130             'm/s             ', 'm/s             ', 'm/s             ', &
1131             'factor          ', 'meters          ', 'meters          ', &
1132             'meters          ', 'none            ', 'none            ', &
1133             'none            ', 'none            ', 'ratio           ', &
1134             'none            ', 'none            ', 'not_used        ', &
1135             'not_used        ', 'not_used        ', 'not_used        ', &
1136             'not_used        ' /)
1137
1138    INTEGER ::  id_dim_prtnum, id_dim_time_pr, id_dim_time_prt, &
1139                id_dim_time_pts, id_dim_time_sp, id_dim_time_ts, id_dim_x_sp, &
1140                id_dim_y_sp, id_dim_zu_sp, id_dim_zw_sp, id_set_pr, &
1141                id_set_prt, id_set_pts, id_set_sp, id_set_ts, id_var_prtnum, &
1142                id_var_rnop_prt, id_var_time_pr, id_var_time_prt, &
1143                id_var_time_pts, id_var_time_sp, id_var_time_ts, id_var_x_sp, &
1144                id_var_y_sp, id_var_zu_sp, id_var_zw_sp, nc_stat
1145
1146    INTEGER, DIMENSION(0:1) ::  id_dim_time_xy, id_dim_time_xz, &
1147                id_dim_time_yz, id_dim_time_3d, id_dim_x_xy, id_dim_xu_xy, &
1148                id_dim_x_xz, id_dim_xu_xz, id_dim_x_yz, id_dim_xu_yz, &
1149                id_dim_x_3d, id_dim_xu_3d, id_dim_y_xy, id_dim_yv_xy, &
1150                id_dim_y_xz, id_dim_yv_xz, id_dim_y_yz, id_dim_yv_yz, &
1151                id_dim_y_3d, id_dim_yv_3d, id_dim_zu_xy, id_dim_zu1_xy, &
1152                id_dim_zu_xz, id_dim_zu_yz, id_dim_zu_3d, id_dim_zw_xy, &
1153                id_dim_zw_xz, id_dim_zw_yz, id_dim_zw_3d, id_set_xy, &
1154                id_set_xz, id_set_yz, id_set_3d, id_var_ind_x_yz, &
1155                id_var_ind_y_xz, id_var_ind_z_xy, id_var_time_xy, &
1156                id_var_time_xz, id_var_time_yz, id_var_time_3d, id_var_x_xy, &
1157                id_var_xu_xy, id_var_x_xz, id_var_xu_xz, id_var_x_yz, &
1158                id_var_xu_yz, id_var_x_3d, id_var_xu_3d, id_var_y_xy, &
1159                id_var_yv_xy, id_var_y_xz, id_var_yv_xz, id_var_y_yz, &
1160                id_var_yv_yz, id_var_y_3d, id_var_yv_3d, id_var_zusi_xy, &
1161                id_var_zusi_3d, id_var_zu_xy, id_var_zu1_xy, id_var_zu_xz, &
1162                id_var_zu_yz, id_var_zu_3d, id_var_zwwi_xy, id_var_zwwi_3d, &
1163                id_var_zw_xy, id_var_zw_xz, id_var_zw_yz, id_var_zw_3d
1164
1165    INTEGER, DIMENSION(10)  ::  id_var_dospx, id_var_dospy
1166    INTEGER, DIMENSION(20)  ::  id_var_prt
1167    INTEGER, DIMENSION(11)  ::  nc_precision
1168    INTEGER, DIMENSION(dopr_norm_num) ::  id_var_norm_dopr
1169
1170    INTEGER, DIMENSION(dopts_num,0:10) ::  id_var_dopts
1171    INTEGER, DIMENSION(0:1,100)        ::  id_var_do2d, id_var_do3d
1172    INTEGER, DIMENSION(100,0:9)        ::  id_dim_z_pr, id_var_dopr, &
1173                                           id_var_z_pr
1174    INTEGER, DIMENSION(dots_max,0:9)   ::  id_var_dots
1175
1176!
1177!-- masked output
1178    CHARACTER (LEN=7), DIMENSION(max_masks,0:1,100) ::  domask_unit
1179
1180    LOGICAL ::  output_for_t0 = .FALSE.
1181
1182    INTEGER, DIMENSION(1:max_masks,0:1) ::  id_dim_time_mask, id_dim_x_mask, &
1183                id_dim_xu_mask, id_dim_y_mask, id_dim_yv_mask, id_dim_zu_mask, &
1184                id_dim_zw_mask, &
1185                id_set_mask, &
1186                id_var_time_mask, id_var_x_mask, id_var_xu_mask, &
1187                id_var_y_mask, id_var_yv_mask, id_var_zu_mask, id_var_zw_mask, &
1188                id_var_zusi_mask, id_var_zwwi_mask
1189
1190    INTEGER, DIMENSION(1:max_masks,0:1,100)         ::  id_var_domask
1191
1192    SAVE
1193
1194 END MODULE netcdf_control
1195
1196
1197
1198 MODULE particle_attributes
1199
1200!------------------------------------------------------------------------------!
1201! Description:
1202! ------------
1203! Definition of variables used to compute particle transport
1204!------------------------------------------------------------------------------!
1205
1206    USE array_kind
1207
1208    CHARACTER (LEN=15)  ::  bc_par_lr = 'cyclic',  bc_par_ns = 'cyclic', &
1209                            bc_par_b  = 'reflect', bc_par_t  = 'absorb', &
1210                            collision_kernel = 'none'
1211
1212#if defined( __parallel )
1213    INTEGER ::  mpi_particle_type
1214#endif
1215    INTEGER ::  deleted_particles = 0, deleted_tails = 0,                      &
1216                dissipation_classes = 10, ibc_par_lr,                          &
1217                ibc_par_ns, ibc_par_b, ibc_par_t, iran_part = -1234567,        &
1218                maximum_number_of_particles = 1000,                            &
1219                maximum_number_of_tailpoints = 100,                            &
1220                maximum_number_of_tails = 0,                                   &
1221                number_of_initial_particles = 0, number_of_particles = 0,      &
1222                number_of_particle_groups = 1, number_of_tails = 0,            &
1223                number_of_initial_tails = 0, offset_ocean_nzt = 0,             &
1224                offset_ocean_nzt_m1 = 0, particles_per_point = 1,              &
1225                particle_file_count = 0, radius_classes = 20,                  &
1226                skip_particles_for_tail = 100, sort_count = 0,                 &
1227                total_number_of_particles, total_number_of_tails = 0,          &
1228                trlp_count_sum, trlp_count_recv_sum, trrp_count_sum,           &
1229                trrp_count_recv_sum, trsp_count_sum, trsp_count_recv_sum,      &
1230                trnp_count_sum, trnp_count_recv_sum
1231
1232    INTEGER, PARAMETER ::  max_number_of_particle_groups = 10
1233
1234    INTEGER, DIMENSION(:), ALLOCATABLE     ::  new_tail_id
1235    INTEGER, DIMENSION(:,:,:), ALLOCATABLE ::  prt_count, prt_start_index
1236
1237    LOGICAL ::  hall_kernel = .FALSE., palm_kernel = .FALSE.,                  &
1238                particle_advection = .FALSE., random_start_position = .FALSE., &
1239                read_particles_from_restartfile = .TRUE.,                      &
1240                uniform_particles = .TRUE., use_kernel_tables = .FALSE.,       &
1241                use_particle_tails = .FALSE., use_sgs_for_particles = .FALSE., &
1242                wang_kernel = .FALSE., write_particle_statistics = .FALSE.
1243
1244    LOGICAL, DIMENSION(max_number_of_particle_groups) ::                       &
1245                vertical_particle_advection = .TRUE.
1246
1247    LOGICAL, DIMENSION(:), ALLOCATABLE ::  particle_mask, tail_mask
1248
1249    REAL    ::  c_0 = 3.0, dt_min_part = 0.0002, dt_prel = 9999999.9,          &
1250                dt_sort_particles = 0.0, dt_write_particle_data = 9999999.9,   &
1251                dvrp_psize = 9999999.9, end_time_prel = 9999999.9,             &
1252                initial_weighting_factor = 1.0,                                &
1253                maximum_tailpoint_age = 100000.0,                              &
1254                minimum_tailpoint_distance = 0.0,                              &
1255                particle_advection_start = 0.0, sgs_wfu_part = 0.3333333,      &
1256                sgs_wfv_part = 0.3333333, sgs_wfw_part = 0.3333333,            &
1257                time_prel = 0.0, time_sort_particles = 0.0,                    &
1258                time_write_particle_data = 0.0
1259
1260    REAL, DIMENSION(max_number_of_particle_groups) ::  &
1261                density_ratio = 9999999.9, pdx = 9999999.9, pdy = 9999999.9, &
1262                pdz = 9999999.9, psb = 9999999.9, psl = 9999999.9,           &
1263                psn = 9999999.9, psr = 9999999.9, pss = 9999999.9,           &
1264                pst = 9999999.9, radius = 9999999.9
1265
1266    REAL, DIMENSION(:,:,:), ALLOCATABLE ::  particle_tail_coordinates
1267
1268
1269    TYPE particle_type
1270       SEQUENCE
1271       REAL    ::  age, age_m, dt_sum, dvrp_psize, e_m, origin_x, origin_y, &
1272                   origin_z, radius, rvar1, rvar2, rvar3, speed_x, speed_y, &
1273                   speed_z, weight_factor, x, y, z
1274       INTEGER ::  class, group, tailpoints, tail_id
1275    END TYPE particle_type
1276
1277    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  initial_particles
1278    TYPE(particle_type), DIMENSION(:), ALLOCATABLE, TARGET ::  part_1, part_2
1279    TYPE(particle_type), DIMENSION(:), POINTER  ::  particles
1280
1281    TYPE particle_groups_type
1282       SEQUENCE
1283       REAL    ::  density_ratio, radius, exp_arg, exp_term
1284    END TYPE particle_groups_type
1285
1286    TYPE(particle_groups_type), DIMENSION(max_number_of_particle_groups) ::&
1287                   particle_groups
1288
1289    SAVE
1290
1291 END MODULE particle_attributes
1292
1293
1294
1295
1296
1297 MODULE pegrid
1298
1299!------------------------------------------------------------------------------!
1300! Description:
1301! ------------
1302! Definition of variables which define processor topology and the exchange of
1303! ghost point layers. This modules must be placed in all routines which contain
1304! MPI-calls.
1305!------------------------------------------------------------------------------!
1306
1307#if defined( __parallel ) && ! defined ( __check )
1308#if defined( __lc )
1309    USE MPI
1310#else
1311    INCLUDE "mpif.h"
1312#endif
1313#endif
1314    CHARACTER(LEN=5)       ::  myid_char = ''
1315    INTEGER                ::  id_inflow = 0, id_recycling = 0, myid = 0,      &
1316                               target_id, npex = -1, npey = -1, numprocs = 1,  &
1317                               numprocs_previous_run = -1,                     &
1318                               tasks_per_node = -9999, threads_per_task = 1
1319
1320    INTEGER, DIMENSION(:,:), ALLOCATABLE ::  hor_index_bounds, &
1321                                             hor_index_bounds_previous_run
1322
1323    LOGICAL ::  collective_wait = .FALSE.
1324
1325#if defined( __parallel )
1326#if defined( __mpi2 )
1327    CHARACTER (LEN=MPI_MAX_PORT_NAME) ::  port_name
1328#endif
1329
1330    INTEGER ::  comm1dx, comm1dy, comm2d, comm_inter, comm_palm, gathered_size,&
1331                ierr, myidx, myidy, ndim = 2, ngp_a, ngp_o, ngp_xy, ngp_y,     &
1332                pleft, pnorth, pright, psouth,                                 &
1333                sendrecvcount_xy, sendrecvcount_yz, sendrecvcount_zx,          &
1334                sendrecvcount_zyd, sendrecvcount_yxd, subdomain_size,          &
1335                type_x, type_x_int, type_xy, type_y, type_y_int
1336
1337    INTEGER ::  ibuf(12), pcoord(2), pdims(2)
1338
1339#if ! defined ( __check )
1340    INTEGER ::  status(MPI_STATUS_SIZE)
1341#endif
1342
1343
1344    INTEGER, DIMENSION(:), ALLOCATABLE ::  ngp_yz, type_xz, type_yz
1345
1346    LOGICAL ::  left_border_pe  = .FALSE., north_border_pe = .FALSE., &
1347                reorder = .TRUE., right_border_pe = .FALSE.,          &
1348                south_border_pe = .FALSE.
1349
1350    LOGICAL, DIMENSION(2) ::  cyclic = (/ .TRUE. , .TRUE. /), &
1351                              remain_dims
1352#endif
1353
1354    SAVE
1355
1356 END MODULE pegrid
1357
1358
1359
1360
1361 MODULE profil_parameter
1362
1363!------------------------------------------------------------------------------!
1364! Description:
1365! ------------
1366! Definition of variables which control PROFIL-output
1367!------------------------------------------------------------------------------!
1368
1369    INTEGER, PARAMETER ::  crmax = 100
1370
1371    CHARACTER (LEN=10), DIMENSION(100) ::  dopr_label = ' '
1372
1373    CHARACTER (LEN=10), DIMENSION(crmax) ::  cross_normalized_x = ' ', &
1374                                             cross_normalized_y = ' '
1375
1376    CHARACTER (LEN=20), DIMENSION(20) ::  cross_ts_profiles = &
1377                           (/ ' E E*               ', ' dt                 ', &
1378                              ' u* w*              ', ' th*                ', &
1379                              ' umax vmax wmax     ', ' div_old div_new    ', &
1380                              ' z_i_wpt z_i_pt     ', ' w"pt"0 w"pt" wpt   ', &
1381                              ' pt(0) pt(zp)       ', ' splux spluy spluz  ', &
1382                              ' L                  ',                         &
1383                            ( '                    ', i9 = 1, 9 ) /)
1384
1385    CHARACTER (LEN=40), DIMENSION(crmax) ::  cross_xtext = &
1386                           (/ 'windspeed in ms>->1         ', &
1387                              'pot. temperature in K       ', &
1388                              'heat flux in K ms>->1       ', &
1389                              'momentum flux in m>2s>2     ', &
1390                              'eddy diffusivity in m>2s>->1', &
1391                              'mixing length in m          ', &
1392                            ( '                            ', i9 = 1, 94 ) /)
1393
1394    CHARACTER (LEN=100), DIMENSION(crmax) ::  cross_profiles = &
1395                           (/ ' u v                           ', &
1396                              ' pt                            ', &
1397                              ' w"pt" w*pt* w*pt*BC wpt wptBC ', &
1398                              ' w"u" w*u* wu w"v" w*v* wv     ', &
1399                              ' km kh                         ', &
1400                              ' l                             ', &
1401                         ( '                               ', i9 = 1, 94 ) /)
1402
1403    INTEGER ::  profile_columns = 3, profile_rows = 2, profile_number = 0
1404
1405    INTEGER ::  cross_linecolors(100,crmax) = 1, &
1406                cross_linestyles(100,crmax) = 0, &
1407                cross_profile_numbers(100,crmax) = 0, &
1408                cross_pnc_local(crmax), cross_profile_number_count(crmax) = 0, &
1409                cross_ts_numbers(crmax,crmax) = 0, &
1410                cross_ts_number_count(crmax) = 0, dopr_crossindex(100) = 0, &
1411                dopr_index(300) = 0, dopr_initial_index(300) = 0, &
1412                dots_crossindex(100) = 0, dots_index(100) = 0, &
1413                linecolors(10) = (/ 2, 3, 4,  5, 7, 8, 12, 15, 16, 23 /), &
1414                linestyles(11) = (/ 0, 7, 3, 10, 4, 1,  9,  2,  5,  8, 6 /)
1415               
1416
1417    REAL ::  cross_normx_factor(100,crmax) = 1.0, &
1418             cross_normy_factor(100,crmax) = 1.0, &
1419             cross_ts_uymax(20) = &
1420                             (/ 999.999, 999.999, 999.999, 999.999, 999.999,   &
1421                                999.999, 999.999, 999.999, 999.999, 999.999,   &
1422                                999.999, 999.999, 999.999, 999.999, 999.999,   &
1423                                999.999, 999.999, 999.999, 999.999, 999.999 /),&
1424             cross_ts_uymax_computed(20) = 999.999, &
1425             cross_ts_uymin(20) = &
1426                             (/ 999.999, 999.999, 999.999,  -5.000, 999.999,   &
1427                                999.999,   0.000, 999.999, 999.999, 999.999,   &
1428                                999.999, 999.999, 999.999, 999.999, 999.999,   &
1429                                999.999, 999.999, 999.999, 999.999, 999.999 /),&
1430             cross_ts_uymin_computed(20) = 999.999, &
1431             cross_uxmax(crmax) = 0.0, cross_uxmax_computed(crmax) = -1.0, &
1432             cross_uxmax_normalized(crmax) = 0.0, &
1433             cross_uxmax_normalized_computed(crmax) = -1.0, &
1434             cross_uxmin(crmax) = 0.0, cross_uxmin_computed(crmax) =  1.0, &
1435             cross_uxmin_normalized(crmax) = 0.0, &
1436             cross_uxmin_normalized_computed(crmax) = 1.0, &
1437             cross_uymax(crmax), cross_uymin(crmax)
1438
1439    SAVE
1440
1441 END MODULE profil_parameter
1442
1443
1444
1445
1446 MODULE spectrum
1447
1448!------------------------------------------------------------------------------!
1449! Description:
1450! ------------
1451! Definition of quantities used for computing spectra
1452!------------------------------------------------------------------------------!
1453
1454    CHARACTER (LEN=6),  DIMENSION(1:5) ::  header_char = (/ 'PS(u) ', 'PS(v) ',&
1455                                           'PS(w) ', 'PS(pt)', 'PS(q) ' /)
1456    CHARACTER (LEN=2),  DIMENSION(10)  ::  spectra_direction = 'x'
1457    CHARACTER (LEN=10), DIMENSION(10)  ::  data_output_sp  = ' '
1458    CHARACTER (LEN=25), DIMENSION(1:5) ::  utext_char =                    &
1459                                           (/ '-power spectrum of u     ', &
1460                                              '-power spectrum of v     ', &
1461                                              '-power spectrum of w     ', &
1462                                              '-power spectrum of ^1185 ', &
1463                                              '-power spectrum of q     ' /)
1464    CHARACTER (LEN=39), DIMENSION(1:5) ::  ytext_char =                        &
1465                                 (/ 'k ^2236 ^2566^2569<u(k) in m>2s>->2    ', &
1466                                    'k ^2236 ^2566^2569<v(k) in m>2s>->2    ', &
1467                                    'k ^2236 ^2566^2569<w(k) in m>2s>->2    ', &
1468                                    'k ^2236 ^2566^2569<^1185(k) in m>2s>->2', &
1469                                    'k ^2236 ^2566^2569<q(k) in m>2s>->2    ' /)
1470
1471    INTEGER ::  klist_x = 0, klist_y = 0, n_sp_x = 0, n_sp_y = 0
1472
1473    INTEGER ::  comp_spectra_level(100) = 999999,                   &
1474                lstyles(100) = (/ 0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1475                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1476                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1477                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1478                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1479                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1480                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1481                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1482                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
1483                                  0, 7, 3, 10, 1, 4, 9, 2, 6, 8 /), &
1484                plot_spectra_level(100) = 999999
1485
1486    REAL    ::  time_to_start_sp = 0.0
1487
1488    SAVE
1489
1490 END MODULE spectrum
1491
1492
1493
1494
1495 MODULE statistics
1496
1497!------------------------------------------------------------------------------!
1498! Description:
1499! ------------
1500! Definition of statistical quantities, e.g. global sums
1501!------------------------------------------------------------------------------!
1502
1503    CHARACTER (LEN=40) ::  region(0:9)
1504    INTEGER ::  pr_palm = 80, statistic_regions = 0
1505    INTEGER ::  u_max_ijk(3), v_max_ijk(3), w_max_ijk(3)
1506    LOGICAL ::  flow_statistics_called = .FALSE.
1507    REAL ::     u_max, v_max, w_max
1508    REAL, DIMENSION(:), ALLOCATABLE       ::  sums_divnew_l, sums_divold_l,   &
1509                                              weight_substep, weight_pres
1510    REAL, DIMENSION(:,:), ALLOCATABLE     ::  sums, sums_wsts_bc_l, ts_value, &
1511                                              sums_wsus_ws_l, sums_wsvs_ws_l, &
1512                                              sums_us2_ws_l, sums_vs2_ws_l,   &
1513                                              sums_ws2_ws_l,                  &
1514                                              sums_wspts_ws_l,                &
1515                                              sums_wssas_ws_l,sums_wsqs_ws_l
1516    REAL, DIMENSION(:,:,:), ALLOCATABLE   ::  hom_sum, rmask, spectrum_x, &
1517                                              spectrum_y, sums_l, sums_l_l, &
1518                                              sums_up_fraction_l
1519    REAL, DIMENSION(:,:,:,:), ALLOCATABLE ::  hom
1520
1521    SAVE
1522
1523 END MODULE statistics
1524
1525
1526
1527
1528 MODULE transpose_indices
1529
1530!------------------------------------------------------------------------------!
1531! Description:
1532! ------------
1533! Definition of indices for transposed arrays
1534!------------------------------------------------------------------------------!
1535
1536    INTEGER ::  nxl_y, nxl_yd, nxl_z, nxr_y, nxr_ya, nxr_yd, nxr_yda, nxr_z, &
1537                nxr_za, nyn_x, nyn_xa, nyn_z, nyn_za, nys_x, nys_z, nzb_x,   &
1538                nzb_y, nzb_yd, nzt_x, nzt_xa, nzt_y, nzt_ya, nzt_yd, nzt_yda
1539               
1540
1541    SAVE
1542
1543 END MODULE transpose_indices
Note: See TracBrowser for help on using the repository browser.