Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |