Home

Context Navigation

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

Last change on this file since 1332 was 1329, checked in by raasch, 11 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 71.3 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |