Home

Context Navigation

data_output_2d.f90 @ 667

Last change on this file since 667 was 667, checked in by suehring, 14 years ago

summary:

Gryschka:

Coupling with different resolution and different numbers of PEs in ocean and atmosphere is available
Exchange of u and v from ocean surface to atmosphere surface
Mirror boundary condition for u and v at the bottom are replaced by dirichlet boundary conditions
Inflow turbulence is now defined by flucuations around spanwise mean
Bugfixes for cyclic_fill and constant_volume_flow

Suehring:

New advection added ( Wicker and Skamarock 5th order ), therefore:
- New module advec_ws.f90
- Modified exchange of ghost boundaries.
- Modified evaluation of turbulent fluxes
- New index bounds nxlg, nxrg, nysg, nyng

advec_ws.f90

Advection scheme for scalars and momentum using the flux formulation of
Wicker and Skamarock 5th order.
Additionally the module contains of a routine using for initialisation and
steering of the statical evaluation. The computation of turbulent fluxes takes
place inside the advection routines.
In case of vector architectures Dirichlet and Radiation boundary conditions are
outstanding and not available. Furthermore simulations within topography are
not possible so far. A further routine local_diss_ij is available and is used
if a control of dissipative fluxes is desired.

check_parameters.f90

Exchange of parameters between ocean and atmosphere via PE0
Check for illegal combination of ws-scheme and timestep scheme.
Check for topography and ws-scheme.
Check for not cyclic boundary conditions in combination with ws-scheme and
loop_optimization = 'vector'.
Check for call_psolver_at_all_substeps and ws-scheme for momentum_advec.

Different processor/grid topology in atmosphere and ocean is now allowed!
Bugfixes in checking for conserve_volume_flow_mode.

exchange_horiz.f90

Dynamic exchange of ghost points with nbgp_local to ensure that no useless
ghost points exchanged in case of multigrid. type_yz(0) and type_xz(0) used for
normal grid, the remaining types used for the several grid levels.
Exchange is done via MPI-Vectors with a dynamic value of ghost points which
depend on the advection scheme. Exchange of left and right PEs is 10% faster
with MPI-Vectors than without.

flow_statistics.f90

When advection is computed with ws-scheme, turbulent fluxes are already
computed in the respective advection routines and buffered in arrays
sums_xxxx_ws_l(). This is due to a consistent treatment of statistics
with the numerics and to avoid unphysical kinks near the surface. So some if-
requests has to be done to dicern between fluxes from ws-scheme other advection
schemes. Furthermore the computation of z_i is only done if the heat flux
exceeds a minimum value. This affects only simulations of a neutral boundary
layer and is due to reasons of computations in the advection scheme.

inflow_turbulence.f90

Using nbgp recycling planes for a better resolution of the turbulent flow near
the inflow.

init_grid.f90

Definition of new array bounds nxlg, nxrg, nysg, nyng on each PE.
Furthermore the allocation of arrays and steering of loops is done with these
parameters. Call of exchange_horiz are modified.
In case of dirichlet bounday condition at the bottom zu(0)=0.0
dzu_mg has to be set explicitly for a equally spaced grid near bottom.
ddzu_pres added to use a equally spaced grid near bottom.

init_pegrid.f90

Moved determination of target_id's from init_coupling
Determination of parameters needed for coupling (coupling_topology, ngp_a, ngp_o)
with different grid/processor-topology in ocean and atmosphere

Adaption of ngp_xy, ngp_y to a dynamic number of ghost points.
The maximum_grid_level changed from 1 to 0. 0 is the normal grid, 1 to
maximum_grid_level the grids for multigrid, in which 0 and 1 are normal grids.
This distinction is due to reasons of data exchange and performance for the
normal grid and grids in poismg.
The definition of MPI-Vectors adapted to a dynamic numer of ghost points.
New MPI-Vectors for data exchange between left and right boundaries added.
This is due to reasons of performance (10% faster).

ATTENTION: nnz_x undefined problem still has to be solved!!!!!!!!
TEST OUTPUT (TO BE REMOVED) logging mpi2 ierr values

parin.f90

Steering parameter dissipation_control added in inipar.

Makefile

Module advec_ws added.

Modules

Removed u_nzb_p1_for_vfc and v_nzb_p1_for_vfc

For coupling with different resolution in ocean and atmophere:
+nx_a, +nx_o, ny_a, +ny_o, ngp_a, ngp_o, +total_2d_o, +total_2d_a,
+coupling_topology

Buffer arrays for the left sided advective fluxes added in arrays_3d.
+flux_s_u, +flux_s_v, +flux_s_w, +diss_s_u, +diss_s_v, +diss_s_w,
+flux_s_pt, +diss_s_pt, +flux_s_e, +diss_s_e, +flux_s_q, +diss_s_q,
+flux_s_sa, +diss_s_sa
3d arrays for dissipation control added. (only necessary for vector arch.)
+var_x, +var_y, +var_z, +gamma_x, +gamma_y, +gamma_z
Default of momentum_advec and scalar_advec changed to 'ws-scheme' .
+exchange_mg added in control_parameters to steer the data exchange.
Parameters +nbgp, +nxlg, +nxrg, +nysg, +nyng added in indices.
flag array +boundary_flags added in indices to steer the degradation of order
of the advective fluxes when non-cyclic boundaries are used.
MPI-datatypes +type_y, +type_y_int and +type_yz for data_exchange added in
pegrid.
+sums_wsus_ws_l, +sums_wsvs_ws_l, +sums_us2_ws_l, +sums_vs2_ws_l,
+sums_ws2_ws_l, +sums_wspts_ws_l, +sums_wssas_ws_l, +sums_wsqs_ws_l
and +weight_substep added in statistics to steer the statistical evaluation
of turbulent fluxes in the advection routines.
LOGICALS +ws_scheme_sca and +ws_scheme_mom added to get a better performance
in prognostic_equations.
LOGICAL +dissipation_control control added to steer numerical dissipation
in ws-scheme.

Changed length of string run_description_header

pres.f90

New allocation of tend when ws-scheme and multigrid is used. This is due to
reasons of perforance of the data_exchange. The same is done with p after
poismg is called.
nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng when no
multigrid is used. Calls of exchange_horiz are modified.

bugfix: After pressure correction no volume flow correction in case of
non-cyclic boundary conditions
(has to be done only before pressure correction)

Call of SOR routine is referenced with ddzu_pres.

prognostic_equations.f90

Calls of the advection routines with WS5 added.
Calls of ws_statistics added to set the statistical arrays to zero after each
time step.

advec_particles.f90

Declaration of de_dx, de_dy, de_dz adapted to additional ghost points.
Furthermore the calls of exchange_horiz were modified.

asselin_filter.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

average_3d_data.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

boundary_conds.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng
Removed mirror boundary conditions for u and v at the bottom in case of
ibc_uv_b == 0. Instead, dirichelt boundary conditions (u=v=0) are set
in init_3d_model

calc_liquid_water_content.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

calc_spectra.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng for
allocation of tend.

check_open.f90

Output of total array size was adapted to nbgp.

data_output_2d.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng in loops and
allocation of arrays local_2d and total_2d.
Calls of exchange_horiz are modified.

data_output_2d.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng in loops and
allocation of arrays. Calls of exchange_horiz are modified.
Skip-value skip_do_avs changed to a dynamic adaption of ghost points.

data_output_mask.f90

Calls of exchange_horiz are modified.

diffusion_e.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

diffusion_s.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

diffusion_u.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

diffusion_v.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

diffusion_w.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

diffusivities.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

diffusivities.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.
Calls of exchange_horiz are modified.

exchange_horiz_2d.f90

Dynamic exchange of ghost points with nbgp, which depends on the advection
scheme. Exchange between left and right PEs is now done with MPI-vectors.

global_min_max.f90

Adapting of the index arrays, because MINLOC assumes lowerbound
at 1 and not at nbgp.

init_3d_model.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng in loops and
allocation of arrays. Calls of exchange_horiz are modified.
Call ws_init to initialize arrays needed for statistical evaluation and
optimization when ws-scheme is used.
Initial volume flow is now calculated by using the variable hom_sum.
Therefore the correction of initial volume flow for non-flat topography
removed (removed u_nzb_p1_for_vfc and v_nzb_p1_for_vfc)
Changed surface boundary conditions for u and v in case of ibc_uv_b == 0 from
mirror bc to dirichlet boundary conditions (u=v=0), so that k=nzb is
representative for the height z0

Bugfix: type conversion of '1' to 64bit for the MAX function (ngp_3d_inner)

init_coupling.f90

determination of target_id's moved to init_pegrid

init_pt_anomaly.f90

Call of exchange_horiz are modified.

init_rankine.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.
Calls of exchange_horiz are modified.

init_slope.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.

header.f90

Output of advection scheme.

poismg.f90

Calls of exchange_horiz are modified.

prandtl_fluxes.f90

Changed surface boundary conditions for u and v from mirror bc to dirichelt bc,
therefore u(uzb,:,:) and v(nzb,:,:) is now representative for the height z0
nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

production_e.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng

read_3d_binary.f90

+/- 1 replaced with +/- nbgp when swapping and allocating variables.

sor.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.
Call of exchange_horiz are modified.
bug removed in declaration of ddzw(), nz replaced by nzt+1

subsidence.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.

sum_up_3d_data.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.

surface_coupler.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng in
MPI_SEND() and MPI_RECV.
additional case for nonequivalent processor and grid topopolgy in ocean and
atmosphere added (coupling_topology = 1)

Added exchange of u and v from Ocean to Atmosphere

time_integration.f90

Calls of exchange_horiz are modified.
Adaption to slooping surface.

timestep.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.

user_3d_data_averaging.f90, user_data_output_2d.f90, user_data_output_3d.f90,
user_actions.f90, user_init.f90, user_init_plant_canopy.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.

user_read_restart_data.f90

Allocation with nbgp.

wall_fluxes.f90

nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng.

write_compressed.f90

Array bounds and nx, ny adapted with nbgp.

sor.f90

bug removed in declaration of ddzw(), nz replaced by nzt+1

Property svn:keywords set to Id

File size: 63.4 KB

Line
1	SUBROUTINE data_output_2d( mode, av )
2
3	!------------------------------------------------------------------------------!
4	! Current revisions:
5	! -----------------
6	! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng in loops and
7	! allocation of arrays local_2d and total_2d.
8	! Calls of exchange_horiz are modiefied.
9	!
10	! Former revisions:
11	! -----------------
12	! $Id: data_output_2d.f90 667 2010-12-23 12:06:00Z suehring $
13	!
14	! 622 2010-12-10 08:08:13Z raasch
15	! optional barriers included in order to speed up collective operations
16	!
17	! 493 2010-03-01 08:30:24Z raasch
18	! NetCDF4 support (parallel output)
19	!
20	! 367 2009-08-25 08:35:52Z maronga
21	! simulated_time in NetCDF output replaced by time_since_reference_point.
22	! Output of NetCDF messages with aid of message handling routine.
23	! Bugfix: averaging along z is not allowed for 2d quantities (e.g. u* and z0)
24	! Output of messages replaced by message handling routine.
25	! Output of user defined 2D (XY) arrays at z=nzb+1 is now possible
26	! Bugfix: to_be_resorted => s_av for time-averaged scalars
27	! Calculation of shf* and qsws* added.
28	!
29	! 215 2008-11-18 09:54:31Z raasch
30	! Bugfix: no output of particle concentration and radius unless particles
31	! have been started
32	!
33	! 96 2007-06-04 08:07:41Z raasch
34	! Output of density and salinity
35	!
36	! 75 2007-03-22 09:54:05Z raasch
37	! Output of precipitation amount/rate and roughness length,
38	! 2nd+3rd argument removed from exchange horiz
39	!
40	! RCS Log replace by Id keyword, revision history cleaned up
41	!
42	! Revision 1.5 2006/08/22 13:50:29 raasch
43	! xz and yz cross sections now up to nzt+1
44	!
45	! Revision 1.2 2006/02/23 10:19:22 raasch
46	! Output of time-averaged data, output of averages along x, y, or z,
47	! output of user-defined quantities,
48	! section data are copied from local_pf to local_2d before they are output,
49	! output of particle concentration and mean radius,
50	! Former subroutine plot_2d renamed data_output_2d, pl2d.. renamed do2d..,
51	! anz renamed ngp, ebene renamed section, pl2d_.._anz renamed do2d_.._n
52	!
53	! Revision 1.1 1997/08/11 06:24:09 raasch
54	! Initial revision
55	!
56	!
57	! Description:
58	! ------------
59	! Data output of horizontal cross-sections in NetCDF format or binary format
60	! compatible to old graphic software iso2d.
61	! Attention: The position of the sectional planes is still not always computed
62	! --------- correctly. (zu is used always)!
63	!------------------------------------------------------------------------------!
64
65	USE arrays_3d
66	USE averaging
67	USE cloud_parameters
68	USE control_parameters
69	USE cpulog
70	USE grid_variables
71	USE indices
72	USE interfaces
73	USE netcdf_control
74	USE particle_attributes
75	USE pegrid
76
77	IMPLICIT NONE
78
79	CHARACTER (LEN=2) :: do2d_mode, mode
80	CHARACTER (LEN=4) :: grid
81	CHARACTER (LEN=25) :: section_chr
82	CHARACTER (LEN=50) :: rtext
83	INTEGER :: av, ngp, file_id, i, if, is, iis, j, k, l, layer_xy, n, psi, &
84	s, sender, &
85	ind(4)
86	LOGICAL :: found, resorted, two_d
87	REAL :: mean_r, s_r3, s_r4
88	REAL, DIMENSION(:), ALLOCATABLE :: level_z
89	REAL, DIMENSION(:,:), ALLOCATABLE :: local_2d, local_2d_l
90	REAL, DIMENSION(:,:,:), ALLOCATABLE :: local_pf
91	#if defined( __parallel )
92	REAL, DIMENSION(:,:), ALLOCATABLE :: total_2d
93	#endif
94	REAL, DIMENSION(:,:,:), POINTER :: to_be_resorted
95
96	NAMELIST /LOCAL/ rtext
97
98	CALL cpu_log (log_point(3),'data_output_2d','start')
99
100	!
101	!-- Immediate return, if no output is requested (no respective sections
102	!-- found in parameter data_output)
103	IF ( mode == 'xy' .AND. .NOT. data_output_xy(av) ) RETURN
104	IF ( mode == 'xz' .AND. .NOT. data_output_xz(av) ) RETURN
105	IF ( mode == 'yz' .AND. .NOT. data_output_yz(av) ) RETURN
106
107	two_d = .FALSE. ! local variable to distinguish between output of pure 2D
108	! arrays and cross-sections of 3D arrays.
109
110	!
111	!-- Depending on the orientation of the cross-section, the respective output
112	!-- files have to be opened.
113	SELECT CASE ( mode )
114
115	CASE ( 'xy' )
116	s = 1
117	ALLOCATE( level_z(nzb:nzt+1), local_2d(nxlg:nxrg,nysg:nyng) )
118
119	!
120	!-- Classic and 64bit offset NetCDF output is done only on PE0.
121	!-- netCDF4/HDF5 output is done in parallel on all PEs.
122	IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
123	THEN
124	CALL check_open( 101+av*10 )
125	ENDIF
126
127	IF ( data_output_2d_on_each_pe ) THEN
128	CALL check_open( 21 )
129	ELSE
130	IF ( myid == 0 ) THEN
131	IF ( iso2d_output ) CALL check_open( 21 )
132	#if defined( __parallel )
133	ALLOCATE( total_2d(-nbgp:nx+nbgp,-nbgp:ny+nbgp) )
134	#endif
135	ENDIF
136	ENDIF
137
138	CASE ( 'xz' )
139	s = 2
140	ALLOCATE( local_2d(nxlg:nxrg,nzb:nzt+1) )
141
142	!
143	!-- Classic and 64bit offset NetCDF output is done only on PE0.
144	!-- netCDF4/HDF5 output may be done in parallel on all PEs.
145	IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
146	THEN
147	CALL check_open( 102+av*10 )
148	ENDIF
149
150	IF ( data_output_2d_on_each_pe ) THEN
151	CALL check_open( 22 )
152	ELSE
153	IF ( myid == 0 ) THEN
154	IF ( iso2d_output ) CALL check_open( 22 )
155	#if defined( __parallel )
156	ALLOCATE( total_2d(-nbgp:nx+nbgp,nzb:nzt+1) )
157	#endif
158	ENDIF
159	ENDIF
160
161	CASE ( 'yz' )
162
163	s = 3
164	ALLOCATE( local_2d(nysg:nyng,nzb:nzt+1) )
165
166	!
167	!-- Classic and 64bit offset NetCDF output is done only on PE0.
168	!-- netCDF4/HDF5 output may be done in parallel on all PEs.
169	IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
170	THEN
171	CALL check_open( 103+av*10 )
172	ENDIF
173
174	IF ( data_output_2d_on_each_pe ) THEN
175	CALL check_open( 23 )
176	ELSE
177	IF ( myid == 0 ) THEN
178	IF ( iso2d_output ) CALL check_open( 23 )
179	#if defined( __parallel )
180	ALLOCATE( total_2d(-nbgp:ny+nbgp,nzb:nzt+1) )
181	#endif
182	ENDIF
183	ENDIF
184
185	CASE DEFAULT
186
187	message_string = 'unknown cross-section: ' // TRIM( mode )
188	CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
189
190	END SELECT
191
192	!
193	!-- Allocate a temporary array for resorting (kji -> ijk).
194	ALLOCATE( local_pf(nxlg:nxrg,nysg:nyng,nzb:nzt+1) )
195
196	!
197	!-- Loop of all variables to be written.
198	!-- Output dimensions chosen
199	if = 1
200	l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
201	do2d_mode = do2d(av,if)(l-1:l)
202
203	DO WHILE ( do2d(av,if)(1:1) /= ' ' )
204
205	IF ( do2d_mode == mode ) THEN
206	!
207	!-- Store the array chosen on the temporary array.
208	resorted = .FALSE.
209	SELECT CASE ( TRIM( do2d(av,if) ) )
210
211	CASE ( 'e_xy', 'e_xz', 'e_yz' )
212	IF ( av == 0 ) THEN
213	to_be_resorted => e
214	ELSE
215	to_be_resorted => e_av
216	ENDIF
217	IF ( mode == 'xy' ) level_z = zu
218
219	CASE ( 'lwp*_xy' ) ! 2d-array
220	IF ( av == 0 ) THEN
221	DO i = nxlg, nxrg
222	DO j = nysg, nyng
223	local_pf(i,j,nzb+1) = SUM( ql(nzb:nzt,j,i) * &
224	dzw(1:nzt+1) )
225	ENDDO
226	ENDDO
227	ELSE
228	DO i = nxlg, nxrg
229	DO j = nysg, nyng
230	local_pf(i,j,nzb+1) = lwp_av(j,i)
231	ENDDO
232	ENDDO
233	ENDIF
234	resorted = .TRUE.
235	two_d = .TRUE.
236	level_z(nzb+1) = zu(nzb+1)
237
238	CASE ( 'p_xy', 'p_xz', 'p_yz' )
239	IF ( av == 0 ) THEN
240	to_be_resorted => p
241	ELSE
242	to_be_resorted => p_av
243	ENDIF
244	IF ( mode == 'xy' ) level_z = zu
245
246	CASE ( 'pc_xy', 'pc_xz', 'pc_yz' ) ! particle concentration
247	IF ( av == 0 ) THEN
248	IF ( simulated_time >= particle_advection_start ) THEN
249	tend = prt_count
250	CALL exchange_horiz( tend, nbgp )
251	ELSE
252	tend = 0.0
253	ENDIF
254	DO i = nxlg, nxrg
255	DO j = nysg, nyng
256	DO k = nzb, nzt+1
257	local_pf(i,j,k) = tend(k,j,i)
258	ENDDO
259	ENDDO
260	ENDDO
261	resorted = .TRUE.
262	ELSE
263	CALL exchange_horiz( pc_av, nbgp )
264	to_be_resorted => pc_av
265	ENDIF
266
267	CASE ( 'pr_xy', 'pr_xz', 'pr_yz' ) ! mean particle radius
268	IF ( av == 0 ) THEN
269	IF ( simulated_time >= particle_advection_start ) THEN
270	DO i = nxl, nxr
271	DO j = nys, nyn
272	DO k = nzb, nzt+1
273	psi = prt_start_index(k,j,i)
274	s_r3 = 0.0
275	s_r4 = 0.0
276	DO n = psi, psi+prt_count(k,j,i)-1
277	s_r3 = s_r3 + particles(n)%radius**3
278	s_r4 = s_r4 + particles(n)%radius**4
279	ENDDO
280	IF ( s_r3 /= 0.0 ) THEN
281	mean_r = s_r4 / s_r3
282	ELSE
283	mean_r = 0.0
284	ENDIF
285	tend(k,j,i) = mean_r
286	ENDDO
287	ENDDO
288	ENDDO
289	CALL exchange_horiz( tend, nbgp )
290	ELSE
291	tend = 0.0
292	END IF
293	DO i = nxlg, nxrg
294	DO j = nysg, nyng
295	DO k = nzb, nzt+1
296	local_pf(i,j,k) = tend(k,j,i)
297	ENDDO
298	ENDDO
299	ENDDO
300	resorted = .TRUE.
301	ELSE
302	CALL exchange_horiz( pr_av, nbgp )
303	to_be_resorted => pr_av
304	ENDIF
305
306	CASE ( 'pra*_xy' ) ! 2d-array / integral quantity => no av
307	CALL exchange_horiz_2d( precipitation_amount )
308	DO i = nxlg, nxrg
309	DO j = nysg, nyng
310	local_pf(i,j,nzb+1) = precipitation_amount(j,i)
311	ENDDO
312	ENDDO
313	precipitation_amount = 0.0 ! reset for next integ. interval
314	resorted = .TRUE.
315	two_d = .TRUE.
316	level_z(nzb+1) = zu(nzb+1)
317
318	CASE ( 'prr*_xy' ) ! 2d-array
319	IF ( av == 0 ) THEN
320	CALL exchange_horiz_2d( precipitation_rate )
321	DO i = nxlg, nxrg
322	DO j = nysg, nyng
323	local_pf(i,j,nzb+1) = precipitation_rate(j,i)
324	ENDDO
325	ENDDO
326	ELSE
327	CALL exchange_horiz_2d( precipitation_rate_av )
328	DO i = nxlg, nxrg
329	DO j = nysg, nyng
330	local_pf(i,j,nzb+1) = precipitation_rate_av(j,i)
331	ENDDO
332	ENDDO
333	ENDIF
334	resorted = .TRUE.
335	two_d = .TRUE.
336	level_z(nzb+1) = zu(nzb+1)
337
338	CASE ( 'pt_xy', 'pt_xz', 'pt_yz' )
339	IF ( av == 0 ) THEN
340	IF ( .NOT. cloud_physics ) THEN
341	to_be_resorted => pt
342	ELSE
343	DO i = nxlg, nxrg
344	DO j = nysg, nyng
345	DO k = nzb, nzt+1
346	local_pf(i,j,k) = pt(k,j,i) + l_d_cp * &
347	pt_d_t(k) * &
348	ql(k,j,i)
349	ENDDO
350	ENDDO
351	ENDDO
352	resorted = .TRUE.
353	ENDIF
354	ELSE
355	to_be_resorted => pt_av
356	ENDIF
357	IF ( mode == 'xy' ) level_z = zu
358
359	CASE ( 'q_xy', 'q_xz', 'q_yz' )
360	IF ( av == 0 ) THEN
361	to_be_resorted => q
362	ELSE
363	to_be_resorted => q_av
364	ENDIF
365	IF ( mode == 'xy' ) level_z = zu
366
367	CASE ( 'ql_xy', 'ql_xz', 'ql_yz' )
368	IF ( av == 0 ) THEN
369	to_be_resorted => ql
370	ELSE
371	to_be_resorted => ql_av
372	ENDIF
373	IF ( mode == 'xy' ) level_z = zu
374
375	CASE ( 'ql_c_xy', 'ql_c_xz', 'ql_c_yz' )
376	IF ( av == 0 ) THEN
377	to_be_resorted => ql_c
378	ELSE
379	to_be_resorted => ql_c_av
380	ENDIF
381	IF ( mode == 'xy' ) level_z = zu
382
383	CASE ( 'ql_v_xy', 'ql_v_xz', 'ql_v_yz' )
384	IF ( av == 0 ) THEN
385	to_be_resorted => ql_v
386	ELSE
387	to_be_resorted => ql_v_av
388	ENDIF
389	IF ( mode == 'xy' ) level_z = zu
390
391	CASE ( 'ql_vp_xy', 'ql_vp_xz', 'ql_vp_yz' )
392	IF ( av == 0 ) THEN
393	to_be_resorted => ql_vp
394	ELSE
395	to_be_resorted => ql_vp_av
396	ENDIF
397	IF ( mode == 'xy' ) level_z = zu
398
399	CASE ( 'qsws*_xy' ) ! 2d-array
400	IF ( av == 0 ) THEN
401	DO i = nxlg, nxrg
402	DO j = nysg, nyng
403	local_pf(i,j,nzb+1) = qsws(j,i)
404	ENDDO
405	ENDDO
406	ELSE
407	DO i = nxlg, nxrg
408	DO j = nysg, nyng
409	local_pf(i,j,nzb+1) = qsws_av(j,i)
410	ENDDO
411	ENDDO
412	ENDIF
413	resorted = .TRUE.
414	two_d = .TRUE.
415	level_z(nzb+1) = zu(nzb+1)
416
417	CASE ( 'qv_xy', 'qv_xz', 'qv_yz' )
418	IF ( av == 0 ) THEN
419	DO i = nxlg, nxrg
420	DO j = nysg, nyng
421	DO k = nzb, nzt+1
422	local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
423	ENDDO
424	ENDDO
425	ENDDO
426	resorted = .TRUE.
427	ELSE
428	to_be_resorted => qv_av
429	ENDIF
430	IF ( mode == 'xy' ) level_z = zu
431
432	CASE ( 'rho_xy', 'rho_xz', 'rho_yz' )
433	IF ( av == 0 ) THEN
434	to_be_resorted => rho
435	ELSE
436	to_be_resorted => rho_av
437	ENDIF
438
439	CASE ( 's_xy', 's_xz', 's_yz' )
440	IF ( av == 0 ) THEN
441	to_be_resorted => q
442	ELSE
443	to_be_resorted => s_av
444	ENDIF
445
446	CASE ( 'sa_xy', 'sa_xz', 'sa_yz' )
447	IF ( av == 0 ) THEN
448	to_be_resorted => sa
449	ELSE
450	to_be_resorted => sa_av
451	ENDIF
452
453	CASE ( 'shf*_xy' ) ! 2d-array
454	IF ( av == 0 ) THEN
455	DO i = nxlg, nxrg
456	DO j = nysg, nyng
457	local_pf(i,j,nzb+1) = shf(j,i)
458	ENDDO
459	ENDDO
460	ELSE
461	DO i = nxlg, nxrg
462	DO j = nysg, nyng
463	local_pf(i,j,nzb+1) = shf_av(j,i)
464	ENDDO
465	ENDDO
466	ENDIF
467	resorted = .TRUE.
468	two_d = .TRUE.
469	level_z(nzb+1) = zu(nzb+1)
470
471	CASE ( 't*_xy' ) ! 2d-array
472	IF ( av == 0 ) THEN
473	DO i = nxlg, nxrg
474	DO j = nysg, nyng
475	local_pf(i,j,nzb+1) = ts(j,i)
476	ENDDO
477	ENDDO
478	ELSE
479	DO i = nxlg, nxrg
480	DO j = nysg, nyng
481	local_pf(i,j,nzb+1) = ts_av(j,i)
482	ENDDO
483	ENDDO
484	ENDIF
485	resorted = .TRUE.
486	two_d = .TRUE.
487	level_z(nzb+1) = zu(nzb+1)
488
489	CASE ( 'u_xy', 'u_xz', 'u_yz' )
490	IF ( av == 0 ) THEN
491	to_be_resorted => u
492	ELSE
493	to_be_resorted => u_av
494	ENDIF
495	IF ( mode == 'xy' ) level_z = zu
496	!
497	!-- Substitute the values generated by "mirror" boundary condition
498	!-- at the bottom boundary by the real surface values.
499	IF ( do2d(av,if) == 'u_xz' .OR. do2d(av,if) == 'u_yz' ) THEN
500	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0
501	ENDIF
502
503	CASE ( 'u*_xy' ) ! 2d-array
504	IF ( av == 0 ) THEN
505	DO i = nxlg, nxrg
506	DO j = nysg, nyng
507	local_pf(i,j,nzb+1) = us(j,i)
508	ENDDO
509	ENDDO
510	ELSE
511	DO i = nxlg, nxrg
512	DO j = nysg, nyng
513	local_pf(i,j,nzb+1) = us_av(j,i)
514	ENDDO
515	ENDDO
516	ENDIF
517	resorted = .TRUE.
518	two_d = .TRUE.
519	level_z(nzb+1) = zu(nzb+1)
520
521	CASE ( 'v_xy', 'v_xz', 'v_yz' )
522	IF ( av == 0 ) THEN
523	to_be_resorted => v
524	ELSE
525	to_be_resorted => v_av
526	ENDIF
527	IF ( mode == 'xy' ) level_z = zu
528	!
529	!-- Substitute the values generated by "mirror" boundary condition
530	!-- at the bottom boundary by the real surface values.
531	IF ( do2d(av,if) == 'v_xz' .OR. do2d(av,if) == 'v_yz' ) THEN
532	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0
533	ENDIF
534
535	CASE ( 'vpt_xy', 'vpt_xz', 'vpt_yz' )
536	IF ( av == 0 ) THEN
537	to_be_resorted => vpt
538	ELSE
539	to_be_resorted => vpt_av
540	ENDIF
541	IF ( mode == 'xy' ) level_z = zu
542
543	CASE ( 'w_xy', 'w_xz', 'w_yz' )
544	IF ( av == 0 ) THEN
545	to_be_resorted => w
546	ELSE
547	to_be_resorted => w_av
548	ENDIF
549	IF ( mode == 'xy' ) level_z = zw
550
551	CASE ( 'z0*_xy' ) ! 2d-array
552	IF ( av == 0 ) THEN
553	DO i = nxlg, nxrg
554	DO j = nysg, nyng
555	local_pf(i,j,nzb+1) = z0(j,i)
556	ENDDO
557	ENDDO
558	ELSE
559	DO i = nxlg, nxrg
560	DO j = nysg, nyng
561	local_pf(i,j,nzb+1) = z0_av(j,i)
562	ENDDO
563	ENDDO
564	ENDIF
565	resorted = .TRUE.
566	two_d = .TRUE.
567	level_z(nzb+1) = zu(nzb+1)
568
569	CASE DEFAULT
570	!
571	!-- User defined quantity
572	CALL user_data_output_2d( av, do2d(av,if), found, grid, &
573	local_pf, two_d )
574	resorted = .TRUE.
575
576	IF ( grid == 'zu' ) THEN
577	IF ( mode == 'xy' ) level_z = zu
578	ELSEIF ( grid == 'zw' ) THEN
579	IF ( mode == 'xy' ) level_z = zw
580	ELSEIF ( grid == 'zu1' ) THEN
581	IF ( mode == 'xy' ) level_z(nzb+1) = zu(nzb+1)
582	ENDIF
583
584	IF ( .NOT. found ) THEN
585	message_string = 'no output provided for: ' // &
586	TRIM( do2d(av,if) )
587	CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
588	ENDIF
589
590	END SELECT
591
592	!
593	!-- Resort the array to be output, if not done above
594	IF ( .NOT. resorted ) THEN
595	DO i = nxlg, nxrg
596	DO j = nysg, nyng
597	DO k = nzb, nzt+1
598	local_pf(i,j,k) = to_be_resorted(k,j,i)
599	ENDDO
600	ENDDO
601	ENDDO
602	ENDIF
603
604	!
605	!-- Output of the individual cross-sections, depending on the cross-
606	!-- section mode chosen.
607	is = 1
608	loop1: DO WHILE ( section(is,s) /= -9999 .OR. two_d )
609
610	SELECT CASE ( mode )
611
612	CASE ( 'xy' )
613	!
614	!-- Determine the cross section index
615	IF ( two_d ) THEN
616	layer_xy = nzb+1
617	ELSE
618	layer_xy = section(is,s)
619	ENDIF
620
621	!
622	!-- Update the NetCDF xy cross section time axis
623	IF ( myid == 0 .OR. netcdf_data_format > 2 ) THEN
624	IF ( simulated_time /= do2d_xy_last_time(av) ) THEN
625	do2d_xy_time_count(av) = do2d_xy_time_count(av) + 1
626	do2d_xy_last_time(av) = simulated_time
627	IF ( ( .NOT. data_output_2d_on_each_pe .AND. &
628	netcdf_output ) .OR. netcdf_data_format > 2 ) &
629	THEN
630	#if defined( __netcdf )
631	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
632	id_var_time_xy(av), &
633	(/ time_since_reference_point /), &
634	start = (/ do2d_xy_time_count(av) /), &
635	count = (/ 1 /) )
636	CALL handle_netcdf_error( 'data_output_2d', 53 )
637	#endif
638	ENDIF
639	ENDIF
640	ENDIF
641	!
642	!-- If required, carry out averaging along z
643	IF ( section(is,s) == -1 .AND. .NOT. two_d ) THEN
644
645	local_2d = 0.0
646	!
647	!-- Carry out the averaging (all data are on the PE)
648	DO k = nzb, nzt+1
649	DO j = nysg, nyng
650	DO i = nxlg, nxrg
651	local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
652	ENDDO
653	ENDDO
654	ENDDO
655
656	local_2d = local_2d / ( nzt -nzb + 2.0)
657
658	ELSE
659	!
660	!-- Just store the respective section on the local array
661	local_2d = local_pf(:,:,layer_xy)
662
663	ENDIF
664
665	#if defined( __parallel )
666	IF ( netcdf_output .AND. netcdf_data_format > 2 ) THEN
667	!
668	!-- Output in NetCDF4/HDF5 format.
669	!-- Do not output redundant ghost point data except for the
670	!-- boundaries of the total domain.
671	IF ( two_d ) THEN
672	iis = 1
673	ELSE
674	iis = is
675	ENDIF
676
677	#if defined( __netcdf )
678	IF ( nxr == nx .AND. nyn /= ny ) THEN
679	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
680	id_var_do2d(av,if), &
681	local_2d(nxl:nxr+1,nys:nyn), &
682	start = (/ nxl+1, nys+1, iis, &
683	do2d_xy_time_count(av) /), &
684	count = (/ nxr-nxl+2, &
685	nyn-nys+1, 1, 1 /) )
686	ELSEIF ( nxr /= nx .AND. nyn == ny ) THEN
687	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
688	id_var_do2d(av,if), &
689	local_2d(nxl:nxr,nys:nyn+1), &
690	start = (/ nxl+1, nys+1, iis, &
691	do2d_xy_time_count(av) /), &
692	count = (/ nxr-nxl+1, &
693	nyn-nys+2, 1, 1 /) )
694	ELSEIF ( nxr == nx .AND. nyn == ny ) THEN
695	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
696	id_var_do2d(av,if), &
697	local_2d(nxl:nxr+1,nys:nyn+1),&
698	start = (/ nxl+1, nys+1, iis, &
699	do2d_xy_time_count(av) /), &
700	count = (/ nxr-nxl+2, &
701	nyn-nys+2, 1, 1 /) )
702	ELSE
703	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
704	id_var_do2d(av,if), &
705	local_2d(nxl:nxr,nys:nyn), &
706	start = (/ nxl+1, nys+1, iis, &
707	do2d_xy_time_count(av) /), &
708	count = (/ nxr-nxl+1, &
709	nyn-nys+1, 1, 1 /) )
710	ENDIF
711
712	CALL handle_netcdf_error( 'data_output_2d', 55 )
713	#endif
714	ELSE
715
716	IF ( data_output_2d_on_each_pe ) THEN
717	!
718	!-- Output of partial arrays on each PE
719	#if defined( __netcdf )
720	IF ( netcdf_output .AND. myid == 0 ) THEN
721	WRITE ( 21 ) simulated_time, &
722	do2d_xy_time_count(av), av
723	ENDIF
724	#endif
725	WRITE ( 21 ) nxlg, nxrg, nysg, nyng
726	WRITE ( 21 ) local_2d
727
728	ELSE
729	!
730	!-- PE0 receives partial arrays from all processors and
731	!-- then outputs them. Here a barrier has to be set,
732	!-- because otherwise "-MPI- FATAL: Remote protocol queue
733	!-- full" may occur.
734	CALL MPI_BARRIER( comm2d, ierr )
735
736	ngp = ( nxrg-nxlg+1 ) * ( nyng-nysg+1 )
737	IF ( myid == 0 ) THEN
738	!
739	!-- Local array can be relocated directly.
740	total_2d(nxlg:nxrg,nysg:nyng) = local_2d
741	!
742	!-- Receive data from all other PEs.
743	DO n = 1, numprocs-1
744	!
745	!-- Receive index limits first, then array.
746	!-- Index limits are received in arbitrary order from
747	!-- the PEs.
748	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
749	MPI_ANY_SOURCE, 0, comm2d, &
750	status, ierr )
751	sender = status(MPI_SOURCE)
752	DEALLOCATE( local_2d )
753	ALLOCATE( local_2d(ind(1):ind(2),ind(3):ind(4)) )
754	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
755	MPI_REAL, sender, 1, comm2d, &
756	status, ierr )
757	total_2d(ind(1):ind(2),ind(3):ind(4)) = local_2d
758	ENDDO
759	!
760	!-- Output of the total cross-section.
761	IF ( iso2d_output ) THEN
762	WRITE (21) total_2d(-nbgp:nx+nbgp,-nbgp:ny+nbgp)
763	ENDIF
764	!
765	!-- Relocate the local array for the next loop increment
766	DEALLOCATE( local_2d )
767	ALLOCATE( local_2d(nxlg:nxrg,nysg:nyng) )
768
769	#if defined( __netcdf )
770	IF ( netcdf_output ) THEN
771	IF ( two_d ) THEN
772	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
773	id_var_do2d(av,if), &
774	total_2d(0:nx+1,0:ny+1), &
775	start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
776	count = (/ nx+2, ny+2, 1, 1 /) )
777	ELSE
778	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
779	id_var_do2d(av,if), &
780	total_2d(0:nx+1,0:ny+1), &
781	start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
782	count = (/ nx+2, ny+2, 1, 1 /) )
783	ENDIF
784	CALL handle_netcdf_error( 'data_output_2d', 54 )
785	ENDIF
786	#endif
787
788	ELSE
789	!
790	!-- First send the local index limits to PE0
791	ind(1) = nxlg; ind(2) = nxrg
792	ind(3) = nysg; ind(4) = nyng
793	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
794	comm2d, ierr )
795	!
796	!-- Send data to PE0
797	CALL MPI_SEND( local_2d(nxlg,nysg), ngp, &
798	MPI_REAL, 0, 1, comm2d, ierr )
799	ENDIF
800	!
801	!-- A barrier has to be set, because otherwise some PEs may
802	!-- proceed too fast so that PE0 may receive wrong data on
803	!-- tag 0
804	CALL MPI_BARRIER( comm2d, ierr )
805	ENDIF
806
807	ENDIF
808	#else
809	IF ( iso2d_output ) THEN
810	WRITE (21) local_2d(nxl:nxr+1,nys:nyn+1)
811	ENDIF
812	#if defined( __netcdf )
813	IF ( netcdf_output ) THEN
814	IF ( two_d ) THEN
815	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
816	id_var_do2d(av,if), &
817	local_2d(nxl:nxr+1,nys:nyn+1), &
818	start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
819	count = (/ nx+2, ny+2, 1, 1 /) )
820	ELSE
821	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
822	id_var_do2d(av,if), &
823	local_2d(nxl:nxr+1,nys:nyn+1), &
824	start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
825	count = (/ nx+2, ny+2, 1, 1 /) )
826	ENDIF
827	CALL handle_netcdf_error( 'data_output_2d', 447 )
828	ENDIF
829	#endif
830	#endif
831	do2d_xy_n = do2d_xy_n + 1
832	!
833	!-- Write LOCAL parameter set for ISO2D.
834	IF ( myid == 0 .AND. iso2d_output ) THEN
835	IF ( section(is,s) /= -1 ) THEN
836	WRITE ( section_chr, '(''z = '',F7.2,'' m (GP '',I3, &
837	&'')'')' &
838	) level_z(layer_xy), layer_xy
839	ELSE
840	section_chr = 'averaged along z'
841	ENDIF
842	IF ( av == 0 ) THEN
843	rtext = TRIM( do2d(av,if) ) // ' t = ' // &
844	TRIM( simulated_time_chr ) // ' ' // &
845	TRIM( section_chr )
846	ELSE
847	rtext = TRIM( do2d(av,if) ) // ' averaged t = ' // &
848	TRIM( simulated_time_chr ) // ' ' // &
849	TRIM( section_chr )
850	ENDIF
851	WRITE (27,LOCAL)
852	ENDIF
853	!
854	!-- For 2D-arrays (e.g. u*) only one cross-section is available.
855	!-- Hence exit loop of output levels.
856	IF ( two_d ) THEN
857	two_d = .FALSE.
858	EXIT loop1
859	ENDIF
860
861	CASE ( 'xz' )
862	!
863	!-- Update the NetCDF xz cross section time axis
864	IF ( myid == 0 .OR. netcdf_data_format > 2 ) THEN
865
866	IF ( simulated_time /= do2d_xz_last_time(av) ) THEN
867	do2d_xz_time_count(av) = do2d_xz_time_count(av) + 1
868	do2d_xz_last_time(av) = simulated_time
869	IF ( ( .NOT. data_output_2d_on_each_pe .AND. &
870	netcdf_output ) .OR. netcdf_data_format > 2 ) &
871	THEN
872	#if defined( __netcdf )
873	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
874	id_var_time_xz(av), &
875	(/ time_since_reference_point /), &
876	start = (/ do2d_xz_time_count(av) /), &
877	count = (/ 1 /) )
878	CALL handle_netcdf_error( 'data_output_2d', 56 )
879	#endif
880	ENDIF
881	ENDIF
882
883	ENDIF
884
885	!
886	!-- If required, carry out averaging along y
887	IF ( section(is,s) == -1 ) THEN
888
889	ALLOCATE( local_2d_l(nxlg:nxrg,nzb:nzt+1) )
890	local_2d_l = 0.0
891	ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
892	!
893	!-- First local averaging on the PE
894	DO k = nzb, nzt+1
895	DO j = nys, nyn
896	DO i = nxlg, nxrg
897	local_2d_l(i,k) = local_2d_l(i,k) + &
898	local_pf(i,j,k)
899	ENDDO
900	ENDDO
901	ENDDO
902	#if defined( __parallel )
903	!
904	!-- Now do the averaging over all PEs along y
905	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
906	CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb), &
907	local_2d(nxlg,nzb), ngp, MPI_REAL, &
908	MPI_SUM, comm1dy, ierr )
909	#else
910	local_2d = local_2d_l
911	#endif
912	local_2d = local_2d / ( ny + 1.0 )
913
914	DEALLOCATE( local_2d_l )
915
916	ELSE
917	!
918	!-- Just store the respective section on the local array
919	!-- (but only if it is available on this PE!)
920	IF ( section(is,s) >= nys .AND. section(is,s) <= nyn ) &
921	THEN
922	local_2d = local_pf(:,section(is,s),nzb:nzt+1)
923	ENDIF
924
925	ENDIF
926
927	#if defined( __parallel )
928	IF ( netcdf_output .AND. netcdf_data_format > 2 ) THEN
929	!
930	!-- ATTENTION: The following lines are a workaround, because
931	!-- independet output does not work with the
932	!-- current NetCDF4 installation. Therefore, data
933	!-- are transferred from PEs having the cross
934	!-- sections to other PEs along y having no cross
935	!-- section data. Some of these data are the
936	!-- output.
937	!-- BEGIN WORKAROUND---------------------------------------
938	IF ( npey /= 1 .AND. section(is,s) /= -1) THEN
939	ALLOCATE( local_2d_l(nxlg:nxrg,nzb:nzt+1) )
940	local_2d_l = 0.0
941	IF ( section(is,s) >= nys .AND. section(is,s) <= nyn )&
942	THEN
943	local_2d_l = local_2d
944	ENDIF
945	#if defined( __parallel )
946	!
947	!-- Distribute data over all PEs along y
948	ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
949	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
950	CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb), &
951	local_2d(nxlg,nzb), ngp, &
952	MPI_REAL, MPI_SUM, comm1dy, ierr )
953	#else
954	local_2d = local_2d_l
955	#endif
956	DEALLOCATE( local_2d_l )
957	ENDIF
958	!-- END WORKAROUND-----------------------------------------
959
960	!
961	!-- Output in NetCDF4/HDF5 format.
962	!-- Output only on those PEs where the respective cross
963	!-- sections reside. Cross sections averaged along y are
964	!-- output on the respective first PE along y (myidy=0).
965	IF ( ( section(is,s) >= nys .AND. &
966	section(is,s) <= nyn ) .OR. &
967	( section(is,s) == -1 .AND. myidy == 0 ) ) THEN
968	!
969	!-- Do not output redundant ghost point data except for the
970	!-- boundaries of the total domain.
971	#if defined( __netcdf )
972	IF ( nxr == nx ) THEN
973	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
974	id_var_do2d(av,if), &
975	local_2d(nxl:nxr+1,nzb:nzt+1), &
976	start = (/ nxl+1, is, 1, &
977	do2d_xz_time_count(av) /), &
978	count = (/ nxr-nxl+2, 1, &
979	nzt+2, 1 /) )
980	ELSE
981	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
982	id_var_do2d(av,if), &
983	local_2d(nxl:nxr,nzb:nzt+1), &
984	start = (/ nxl+1, is, 1, &
985	do2d_xz_time_count(av) /), &
986	count = (/ nxr-nxl+1, 1, &
987	nzt+2, 1 /) )
988	ENDIF
989
990	CALL handle_netcdf_error( 'data_output_2d', 57 )
991
992	ELSE
993	!
994	!-- Output on other PEs. Only one point is output!!
995	!-- ATTENTION: This is a workaround (see above)!!
996	IF ( npey /= 1 ) THEN
997	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
998	id_var_do2d(av,if), &
999	local_2d(nxl:nxl,nzb:nzb), &
1000	start = (/ nxl+1, is, 1, &
1001	do2d_xz_time_count(av) /), &
1002	count = (/ 1, 1, 1, 1 /) )
1003	CALL handle_netcdf_error( 'data_output_2d', 451 )
1004	ENDIF
1005	#endif
1006	ENDIF
1007
1008	ELSE
1009
1010	IF ( data_output_2d_on_each_pe ) THEN
1011	!
1012	!-- Output of partial arrays on each PE. If the cross
1013	!-- section does not reside on the PE, output special
1014	!-- index values.
1015	#if defined( __netcdf )
1016	IF ( netcdf_output .AND. myid == 0 ) THEN
1017	WRITE ( 22 ) simulated_time, &
1018	do2d_xz_time_count(av), av
1019	ENDIF
1020	#endif
1021	IF ( ( section(is,s) >= nys .AND. &
1022	section(is,s) <= nyn ) .OR. &
1023	( section(is,s) == -1 .AND. nys-1 == -1 ) ) &
1024	THEN
1025	WRITE (22) nxlg, nxrg, nzb, nzt+1
1026	WRITE (22) local_2d
1027	ELSE
1028	WRITE (22) -1, -1, -1, -1
1029	ENDIF
1030
1031	ELSE
1032	!
1033	!-- PE0 receives partial arrays from all processors of the
1034	!-- respective cross section and outputs them. Here a
1035	!-- barrier has to be set, because otherwise
1036	!-- "-MPI- FATAL: Remote protocol queue full" may occur.
1037	CALL MPI_BARRIER( comm2d, ierr )
1038
1039	ngp = ( nxrg-nxlg+1 ) * ( nzt-nzb+2 )
1040	IF ( myid == 0 ) THEN
1041	!
1042	!-- Local array can be relocated directly.
1043	IF ( ( section(is,s) >= nys .AND. &
1044	section(is,s) <= nyn ) .OR. &
1045	( section(is,s) == -1 .AND. nys-1 == -1 ) ) &
1046	THEN
1047	total_2d(nxlg:nxrg,nzb:nzt+1) = local_2d
1048	ENDIF
1049	!
1050	!-- Receive data from all other PEs.
1051	DO n = 1, numprocs-1
1052	!
1053	!-- Receive index limits first, then array.
1054	!-- Index limits are received in arbitrary order from
1055	!-- the PEs.
1056	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
1057	MPI_ANY_SOURCE, 0, comm2d, &
1058	status, ierr )
1059	!
1060	!-- Not all PEs have data for XZ-cross-section.
1061	IF ( ind(1) /= -9999 ) THEN
1062	sender = status(MPI_SOURCE)
1063	DEALLOCATE( local_2d )
1064	ALLOCATE( local_2d(ind(1):ind(2), &
1065	ind(3):ind(4)) )
1066	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1067	MPI_REAL, sender, 1, comm2d, &
1068	status, ierr )
1069	total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1070	local_2d
1071	ENDIF
1072	ENDDO
1073	!
1074	!-- Output of the total cross-section.
1075	IF ( iso2d_output ) THEN
1076	WRITE (22) total_2d(-nbgp:nx+nbgp,nzb:nzt+1)
1077	ENDIF
1078	!
1079	!-- Relocate the local array for the next loop increment
1080	DEALLOCATE( local_2d )
1081	ALLOCATE( local_2d(nxlg:nxrg,nzb:nzt+1) )
1082
1083	#if defined( __netcdf )
1084	IF ( netcdf_output ) THEN
1085	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
1086	id_var_do2d(av,if), &
1087	total_2d(0:nx+1,nzb:nzt+1),&
1088	start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1089	count = (/ nx+2, 1, nz+2, 1 /) )
1090	CALL handle_netcdf_error( 'data_output_2d', 58 )
1091	ENDIF
1092	#endif
1093
1094	ELSE
1095	!
1096	!-- If the cross section resides on the PE, send the
1097	!-- local index limits, otherwise send -9999 to PE0.
1098	IF ( ( section(is,s) >= nys .AND. &
1099	section(is,s) <= nyn ) .OR. &
1100	( section(is,s) == -1 .AND. nys-1 == -1 ) ) &
1101	THEN
1102	ind(1) = nxlg; ind(2) = nxrg
1103	ind(3) = nzb; ind(4) = nzt+1
1104	ELSE
1105	ind(1) = -9999; ind(2) = -9999
1106	ind(3) = -9999; ind(4) = -9999
1107	ENDIF
1108	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1109	comm2d, ierr )
1110	!
1111	!-- If applicable, send data to PE0.
1112	IF ( ind(1) /= -9999 ) THEN
1113	CALL MPI_SEND( local_2d(nxlg,nzb), ngp, &
1114	MPI_REAL, 0, 1, comm2d, ierr )
1115	ENDIF
1116	ENDIF
1117	!
1118	!-- A barrier has to be set, because otherwise some PEs may
1119	!-- proceed too fast so that PE0 may receive wrong data on
1120	!-- tag 0
1121	CALL MPI_BARRIER( comm2d, ierr )
1122	ENDIF
1123
1124	ENDIF
1125	#else
1126	IF ( iso2d_output ) THEN
1127	WRITE (22) local_2d(nxl:nxr+1,nzb:nzt+1)
1128	ENDIF
1129	#if defined( __netcdf )
1130	IF ( netcdf_output ) THEN
1131	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
1132	id_var_do2d(av,if), &
1133	local_2d(nxl:nxr+1,nzb:nzt+1), &
1134	start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1135	count = (/ nx+2, 1, nz+2, 1 /) )
1136	CALL handle_netcdf_error( 'data_output_2d', 451 )
1137	ENDIF
1138	#endif
1139	#endif
1140	do2d_xz_n = do2d_xz_n + 1
1141	!
1142	!-- Write LOCAL-parameter set for ISO2D.
1143	IF ( myid == 0 .AND. iso2d_output ) THEN
1144	IF ( section(is,s) /= -1 ) THEN
1145	WRITE ( section_chr, '(''y = '',F8.2,'' m (GP '',I3, &
1146	&'')'')' &
1147	) section(is,s)*dy, section(is,s)
1148	ELSE
1149	section_chr = 'averaged along y'
1150	ENDIF
1151	IF ( av == 0 ) THEN
1152	rtext = TRIM( do2d(av,if) ) // ' t = ' // &
1153	TRIM( simulated_time_chr ) // ' ' // &
1154	TRIM( section_chr )
1155	ELSE
1156	rtext = TRIM( do2d(av,if) ) // ' averaged t = ' // &
1157	TRIM( simulated_time_chr ) // ' ' // &
1158	TRIM( section_chr )
1159	ENDIF
1160	WRITE (28,LOCAL)
1161	ENDIF
1162
1163	CASE ( 'yz' )
1164	!
1165	!-- Update the NetCDF yz cross section time axis
1166	IF ( myid == 0 .OR. netcdf_data_format > 2 ) THEN
1167
1168	IF ( simulated_time /= do2d_yz_last_time(av) ) THEN
1169	do2d_yz_time_count(av) = do2d_yz_time_count(av) + 1
1170	do2d_yz_last_time(av) = simulated_time
1171	IF ( ( .NOT. data_output_2d_on_each_pe .AND. &
1172	netcdf_output ) .OR. netcdf_data_format > 2 ) &
1173	THEN
1174	#if defined( __netcdf )
1175	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1176	id_var_time_yz(av), &
1177	(/ time_since_reference_point /), &
1178	start = (/ do2d_yz_time_count(av) /), &
1179	count = (/ 1 /) )
1180	CALL handle_netcdf_error( 'data_output_2d', 59 )
1181	#endif
1182	ENDIF
1183	ENDIF
1184
1185	ENDIF
1186	!
1187	!-- If required, carry out averaging along x
1188	IF ( section(is,s) == -1 ) THEN
1189
1190	ALLOCATE( local_2d_l(nysg:nyng,nzb:nzt+1) )
1191	local_2d_l = 0.0
1192	ngp = ( nyng-nysg+1 ) * ( nzt-nzb+2 )
1193	!
1194	!-- First local averaging on the PE
1195	DO k = nzb, nzt+1
1196	DO j = nysg, nyng
1197	DO i = nxl, nxr
1198	local_2d_l(j,k) = local_2d_l(j,k) + &
1199	local_pf(i,j,k)
1200	ENDDO
1201	ENDDO
1202	ENDDO
1203	#if defined( __parallel )
1204	!
1205	!-- Now do the averaging over all PEs along x
1206	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1207	CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb), &
1208	local_2d(nysg,nzb), ngp, MPI_REAL, &
1209	MPI_SUM, comm1dx, ierr )
1210	#else
1211	local_2d = local_2d_l
1212	#endif
1213	local_2d = local_2d / ( nx + 1.0 )
1214
1215	DEALLOCATE( local_2d_l )
1216
1217	ELSE
1218	!
1219	!-- Just store the respective section on the local array
1220	!-- (but only if it is available on this PE!)
1221	IF ( section(is,s) >= nxl .AND. section(is,s) <= nxr ) &
1222	THEN
1223	local_2d = local_pf(section(is,s),:,nzb:nzt+1)
1224	ENDIF
1225
1226	ENDIF
1227
1228	#if defined( __parallel )
1229	IF ( netcdf_output .AND. netcdf_data_format > 2 ) THEN
1230	!
1231	!-- ATTENTION: The following lines are a workaround, because
1232	!-- independet output does not work with the
1233	!-- current NetCDF4 installation. Therefore, data
1234	!-- are transferred from PEs having the cross
1235	!-- sections to other PEs along y having no cross
1236	!-- section data. Some of these data are the
1237	!-- output.
1238	!-- BEGIN WORKAROUND---------------------------------------
1239	IF ( npex /= 1 .AND. section(is,s) /= -1) THEN
1240	ALLOCATE( local_2d_l(nysg:nyng,nzb:nzt+1) )
1241	local_2d_l = 0.0
1242	IF ( section(is,s) >= nxl .AND. section(is,s) <= nxr )&
1243	THEN
1244	local_2d_l = local_2d
1245	ENDIF
1246	#if defined( __parallel )
1247	!
1248	!-- Distribute data over all PEs along x
1249	ngp = ( nyng-nysg+1 ) * ( nzt-nzb + 2 )
1250	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1251	CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb), &
1252	local_2d(nysg,nzb), ngp, &
1253	MPI_REAL, MPI_SUM, comm1dx, ierr )
1254	#else
1255	local_2d = local_2d_l
1256	#endif
1257	DEALLOCATE( local_2d_l )
1258	ENDIF
1259	!-- END WORKAROUND-----------------------------------------
1260
1261	!
1262	!-- Output in NetCDF4/HDF5 format.
1263	!-- Output only on those PEs where the respective cross
1264	!-- sections reside. Cross sections averaged along x are
1265	!-- output on the respective first PE along x (myidx=0).
1266	IF ( ( section(is,s) >= nxl .AND. &
1267	section(is,s) <= nxr ) .OR. &
1268	( section(is,s) == -1 .AND. myidx == 0 ) ) THEN
1269	!
1270	!-- Do not output redundant ghost point data except for the
1271	!-- boundaries of the total domain.
1272	#if defined( __netcdf )
1273	IF ( nyn == ny ) THEN
1274	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1275	id_var_do2d(av,if), &
1276	local_2d(nys:nyn+1,nzb:nzt+1), &
1277	start = (/ is, nys+1, 1, &
1278	do2d_yz_time_count(av) /), &
1279	count = (/ 1, nyn-nys+2, &
1280	nzt+2, 1 /) )
1281	ELSE
1282	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1283	id_var_do2d(av,if), &
1284	local_2d(nys:nyn,nzb:nzt+1), &
1285	start = (/ is, nys+1, 1, &
1286	do2d_yz_time_count(av) /), &
1287	count = (/ 1, nyn-nys+1, &
1288	nzt+2, 1 /) )
1289	ENDIF
1290
1291	CALL handle_netcdf_error( 'data_output_2d', 60 )
1292
1293	ELSE
1294	!
1295	!-- Output on other PEs. Only one point is output!!
1296	!-- ATTENTION: This is a workaround (see above)!!
1297	IF ( npex /= 1 ) THEN
1298	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1299	id_var_do2d(av,if), &
1300	local_2d(nys:nys,nzb:nzb), &
1301	start = (/ is, nys+1, 1, &
1302	do2d_yz_time_count(av) /), &
1303	count = (/ 1, 1, 1, 1 /) )
1304	CALL handle_netcdf_error( 'data_output_2d', 452 )
1305	ENDIF
1306	#endif
1307	ENDIF
1308
1309	ELSE
1310
1311	IF ( data_output_2d_on_each_pe ) THEN
1312	!
1313	!-- Output of partial arrays on each PE. If the cross
1314	!-- section does not reside on the PE, output special
1315	!-- index values.
1316	#if defined( __netcdf )
1317	IF ( netcdf_output .AND. myid == 0 ) THEN
1318	WRITE ( 23 ) simulated_time, &
1319	do2d_yz_time_count(av), av
1320	ENDIF
1321	#endif
1322	IF ( ( section(is,s) >= nxl .AND. &
1323	section(is,s) <= nxr ) .OR. &
1324	( section(is,s) == -1 .AND. nxl-1 == -1 ) ) &
1325	THEN
1326	WRITE (23) nysg, nyng, nzb, nzt+1
1327	WRITE (23) local_2d
1328	ELSE
1329	WRITE (23) -1, -1, -1, -1
1330	ENDIF
1331
1332	ELSE
1333	!
1334	!-- PE0 receives partial arrays from all processors of the
1335	!-- respective cross section and outputs them. Here a
1336	!-- barrier has to be set, because otherwise
1337	!-- "-MPI- FATAL: Remote protocol queue full" may occur.
1338	CALL MPI_BARRIER( comm2d, ierr )
1339
1340	ngp = ( nyng-nysg+1 ) * ( nzt-nzb+2 )
1341	IF ( myid == 0 ) THEN
1342	!
1343	!-- Local array can be relocated directly.
1344	IF ( ( section(is,s) >= nxl .AND. &
1345	section(is,s) <= nxr ) .OR. &
1346	( section(is,s) == -1 .AND. nxl-1 == -1 ) ) &
1347	THEN
1348	total_2d(nysg:nyng,nzb:nzt+1) = local_2d
1349	ENDIF
1350	!
1351	!-- Receive data from all other PEs.
1352	DO n = 1, numprocs-1
1353	!
1354	!-- Receive index limits first, then array.
1355	!-- Index limits are received in arbitrary order from
1356	!-- the PEs.
1357	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
1358	MPI_ANY_SOURCE, 0, comm2d, &
1359	status, ierr )
1360	!
1361	!-- Not all PEs have data for YZ-cross-section.
1362	IF ( ind(1) /= -9999 ) THEN
1363	sender = status(MPI_SOURCE)
1364	DEALLOCATE( local_2d )
1365	ALLOCATE( local_2d(ind(1):ind(2), &
1366	ind(3):ind(4)) )
1367	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1368	MPI_REAL, sender, 1, comm2d, &
1369	status, ierr )
1370	total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1371	local_2d
1372	ENDIF
1373	ENDDO
1374	!
1375	!-- Output of the total cross-section.
1376	IF ( iso2d_output ) THEN
1377	WRITE (23) total_2d(0:ny+1,nzb:nzt+1)
1378	ENDIF
1379	!
1380	!-- Relocate the local array for the next loop increment
1381	DEALLOCATE( local_2d )
1382	ALLOCATE( local_2d(nysg:nyng,nzb:nzt+1) )
1383
1384	#if defined( __netcdf )
1385	IF ( netcdf_output ) THEN
1386	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1387	id_var_do2d(av,if), &
1388	total_2d(0:ny+1,nzb:nzt+1),&
1389	start = (/ is, 1, 1, do2d_yz_time_count(av) /), &
1390	count = (/ 1, ny+2, nz+2, 1 /) )
1391	CALL handle_netcdf_error( 'data_output_2d', 61 )
1392	ENDIF
1393	#endif
1394
1395	ELSE
1396	!
1397	!-- If the cross section resides on the PE, send the
1398	!-- local index limits, otherwise send -9999 to PE0.
1399	IF ( ( section(is,s) >= nxl .AND. &
1400	section(is,s) <= nxr ) .OR. &
1401	( section(is,s) == -1 .AND. nxl-1 == -1 ) ) &
1402	THEN
1403	ind(1) = nysg; ind(2) = nyng
1404	ind(3) = nzb; ind(4) = nzt+1
1405	ELSE
1406	ind(1) = -9999; ind(2) = -9999
1407	ind(3) = -9999; ind(4) = -9999
1408	ENDIF
1409	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1410	comm2d, ierr )
1411	!
1412	!-- If applicable, send data to PE0.
1413	IF ( ind(1) /= -9999 ) THEN
1414	CALL MPI_SEND( local_2d(nysg,nzb), ngp, &
1415	MPI_REAL, 0, 1, comm2d, ierr )
1416	ENDIF
1417	ENDIF
1418	!
1419	!-- A barrier has to be set, because otherwise some PEs may
1420	!-- proceed too fast so that PE0 may receive wrong data on
1421	!-- tag 0
1422	CALL MPI_BARRIER( comm2d, ierr )
1423	ENDIF
1424
1425	ENDIF
1426	#else
1427	IF ( iso2d_output ) THEN
1428	WRITE (23) local_2d(nys:nyn+1,nzb:nzt+1)
1429	ENDIF
1430	#if defined( __netcdf )
1431	IF ( netcdf_output ) THEN
1432	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1433	id_var_do2d(av,if), &
1434	local_2d(nys:nyn+1,nzb:nzt+1), &
1435	start = (/ is, 1, 1, do2d_xz_time_count(av) /), &
1436	count = (/ 1, ny+2, nz+2, 1 /) )
1437	CALL handle_netcdf_error( 'data_output_2d', 452 )
1438	ENDIF
1439	#endif
1440	#endif
1441	do2d_yz_n = do2d_yz_n + 1
1442	!
1443	!-- Write LOCAL-parameter set for ISO2D.
1444	IF ( myid == 0 .AND. iso2d_output ) THEN
1445	IF ( section(is,s) /= -1 ) THEN
1446	WRITE ( section_chr, '(''x = '',F8.2,'' m (GP '',I3, &
1447	&'')'')' &
1448	) section(is,s)*dx, section(is,s)
1449	ELSE
1450	section_chr = 'averaged along x'
1451	ENDIF
1452	IF ( av == 0 ) THEN
1453	rtext = TRIM( do2d(av,if) ) // ' t = ' // &
1454	TRIM( simulated_time_chr ) // ' ' // &
1455	TRIM( section_chr )
1456	ELSE
1457	rtext = TRIM( do2d(av,if) ) // ' averaged t = ' // &
1458	TRIM( simulated_time_chr ) // ' ' // &
1459	TRIM( section_chr )
1460	ENDIF
1461	WRITE (29,LOCAL)
1462	ENDIF
1463
1464	END SELECT
1465
1466	is = is + 1
1467	ENDDO loop1
1468
1469	ENDIF
1470
1471	if = if + 1
1472	l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
1473	do2d_mode = do2d(av,if)(l-1:l)
1474
1475	ENDDO
1476
1477	!
1478	!-- Deallocate temporary arrays.
1479	IF ( ALLOCATED( level_z ) ) DEALLOCATE( level_z )
1480	DEALLOCATE( local_pf, local_2d )
1481	#if defined( __parallel )
1482	IF ( .NOT. data_output_2d_on_each_pe .AND. myid == 0 ) THEN
1483	DEALLOCATE( total_2d )
1484	ENDIF
1485	#endif
1486
1487	!
1488	!-- Close plot output file.
1489	file_id = 20 + s
1490
1491	IF ( data_output_2d_on_each_pe ) THEN
1492	CALL close_file( file_id )
1493	ELSE
1494	IF ( myid == 0 ) CALL close_file( file_id )
1495	ENDIF
1496
1497
1498	CALL cpu_log (log_point(3),'data_output_2d','stop','nobarrier')
1499
1500	END SUBROUTINE data_output_2d

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: palm/trunk/SOURCE/data_output_2d.f90 @ 667

Download in other formats: