Home

Context Navigation

data_output_3d.f90 @ 4867

Last change on this file since 4867 was 4828, checked in by Giersch, 4 years ago
Copyright updated to year 2021, interface pmc_sort removed to accelarate the nesting code
Property svn:keywords set to `Id`
File size: 26.6 KB

Line
1	!> @file data_output_3d.f90
2	!--------------------------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
6	! Public License as published by the Free Software Foundation, either version 3 of the License, or
7	! (at your option) any later version.
8	!
9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
10	! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11	! Public License for more details.
12	!
13	! You should have received a copy of the GNU General Public License along with PALM. If not, see
14	! <http://www.gnu.org/licenses/>.
15	!
16	! Copyright 1997-2021 Leibniz Universitaet Hannover
17	!--------------------------------------------------------------------------------------------------!
18	!
19	! Current revisions:
20	! ------------------
21	!
22	!
23	! Former revisions:
24	! -----------------
25	! $Id: data_output_3d.f90 4828 2021-01-05 11:21:41Z moh.hefny $
26	! Keep 32-bit output for netcdf_data_output < 5
27	!
28	! 4768 2020-11-02 19:11:23Z suehring
29	! Enable 3D data output also with 64-bit precision
30	!
31	! 4559 2020-06-11 08:51:48Z raasch
32	! file re-formatted to follow the PALM coding standard
33	!
34	! 4457 2020-03-11 14:20:43Z raasch
35	! use statement for exchange horiz added
36	!
37	! 4444 2020-03-05 15:59:50Z raasch
38	! bugfix: cpp-directives for serial mode added
39	!
40	! 4360 2020-01-07 11:25:50Z suehring
41	! Introduction of wall_flags_total_0, which currently sets bits based on static topography
42	! information used in wall_flags_static_0
43	!
44	! 4329 2019-12-10 15:46:36Z motisi
45	! Renamed wall_flags_0 to wall_flags_static_0
46	!
47	! 4182 2019-08-22 15:20:23Z scharf
48	! Corrected "Former revisions" section
49	!
50	! 4162 2019-08-16 05:54:29Z raasch
51	! bugfix for r4155
52	!
53	! 4155 2019-08-14 06:25:18Z raasch
54	! bugfix for 3d-output in serial mode (ghost points must not be written)
55	!
56	! 4127 2019-07-30 14:47:10Z suehring
57	! Adjustment for top boundary index for plant-canopy model outputs (merge from branch resler)
58	!
59	! 4048 2019-06-21 21:00:21Z knoop
60	! Moved tcm_data_output_3d to module_interface
61	!
62	! 4039 2019-06-18 10:32:41Z suehring
63	! modularize diagnostic output
64	!
65	! 3994 2019-05-22 18:08:09Z suehring
66	! output of turbulence intensity added
67	!
68	! 3987 2019-05-22 09:52:13Z kanani
69	! Introduce alternative switch for debug output during timestepping
70	!
71	! 3885 2019-04-11 11:29:34Z kanani
72	! Changes related to global restructuring of location messages and introduction of additional debug
73	! messages
74	!
75	! 3814 2019-03-26 08:40:31Z pavelkrc
76	! unused variables removed
77	!
78	! 3655 2019-01-07 16:51:22Z knoop
79	! Bugfix: use time_since_reference_point instead of simulated_time (relevant when using wall/soil
80	! spinup)
81	!
82	! Revision 1.1 1997/09/03 06:29:36 raasch
83	! Initial revision
84	!
85	!
86	! Description:
87	! ------------
88	!> Output of the 3D-arrays in netCDF and/or AVS format.
89	!--------------------------------------------------------------------------------------------------!
90	SUBROUTINE data_output_3d( av )
91
92
93	USE arrays_3d, &
94	ONLY: d_exner, e, p, pt, q, ql, ql_c, ql_v, s, tend, u, v, vpt, w
95
96	USE averaging
97
98	USE basic_constants_and_equations_mod, &
99	ONLY: lv_d_cp
100
101	USE bulk_cloud_model_mod, &
102	ONLY: bulk_cloud_model
103
104	USE control_parameters, &
105	ONLY: debug_output_timestep, do3d, do3d_no, do3d_time_count, land_surface, message_string,&
106	ntdim_3d, nz_do3d, plant_canopy, psolver, time_since_reference_point, urban_surface,&
107	varnamelength
108
109	#if defined( __parallel )
110	USE control_parameters, &
111	ONLY: io_blocks, io_group
112	#endif
113
114	USE cpulog, &
115	ONLY: cpu_log, log_point
116
117	USE exchange_horiz_mod, &
118	ONLY: exchange_horiz
119
120	USE indices, &
121	ONLY: nbgp, nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt, wall_flags_total_0
122
123	#if ! defined( __parallel )
124	USE indices, &
125	ONLY: nx, ny
126	#endif
127
128	USE kinds
129
130	USE land_surface_model_mod, &
131	ONLY: lsm_data_output_3d, nzb_soil, nzt_soil
132
133	USE module_interface, &
134	ONLY: module_interface_data_output_3d
135
136	#if defined( __netcdf )
137	USE NETCDF
138	#endif
139
140	USE netcdf_interface, &
141	ONLY: fill_value, id_set_3d, id_var_do3d, id_var_time_3d, nc_stat, netcdf_data_format, &
142	netcdf_handle_error
143
144	USE particle_attributes, &
145	ONLY: grid_particles, number_of_particles, particles, particle_advection_start, prt_count
146
147	USE pegrid
148
149	USE plant_canopy_model_mod, &
150	ONLY: pch_index
151
152	USE radiation_model_mod, &
153	ONLY: nz_urban_b, nz_urban_t
154
155	USE urban_surface_mod, &
156	ONLY: usm_data_output_3d
157
158
159	IMPLICIT NONE
160
161	INTEGER(iwp) :: av !< flag for (non-)average output
162	INTEGER(iwp) :: flag_nr !< number of masking flag
163	INTEGER(iwp) :: i !< loop index
164	INTEGER(iwp) :: ivar !< variable index
165	INTEGER(iwp) :: j !< loop index
166	INTEGER(iwp) :: k !< loop index
167	INTEGER(iwp) :: n !< loop index
168	INTEGER(iwp) :: nzb_do !< vertical lower limit for data output
169	INTEGER(iwp) :: nzt_do !< vertical upper limit for data output
170
171	LOGICAL :: found !< true if output variable was found
172	LOGICAL :: resorted !< true if variable is resorted
173
174	REAL(wp) :: mean_r !< mean particle radius
175	REAL(wp) :: s_r2 !< sum( particle-radius**2 )
176	REAL(wp) :: s_r3 !< sum( particle-radius**3 )
177
178	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: local_pf !< output array
179
180	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< pointer to array which shall be
181	!< output
182
183	CHARACTER (LEN=varnamelength) :: trimvar !< TRIM of output-variable string
184
185	!
186	!-- Return, if nothing to output
187	IF ( do3d_no(av) == 0 ) RETURN
188
189	IF ( debug_output_timestep ) CALL debug_message( 'data_output_3d', 'start' )
190
191	CALL cpu_log (log_point(14),'data_output_3d','start')
192
193	!
194	!-- Open output file.
195	!-- For classic or 64bit netCDF output on more than one PE, each PE opens its own file and writes
196	!-- the data of its subdomain in binary format. After the run, these files are combined to one
197	!-- NetCDF file by combine_plot_fields.
198	!-- For netCDF4/HDF5 output, data is written in parallel into one file.
199	IF ( netcdf_data_format < 5 ) THEN
200	#if defined( __parallel )
201	CALL check_open( 30 )
202	#endif
203	IF ( myid == 0 ) CALL check_open( 106+av*10 )
204	ELSE
205	CALL check_open( 106+av*10 )
206	ENDIF
207
208	!
209	!-- For parallel netcdf output the time axis must be limited. Return, if this limit is exceeded.
210	!-- This could be the case, if the simulated time exceeds the given end time by the length of the
211	!-- given output interval.
212	IF ( netcdf_data_format > 4 ) THEN
213	IF ( do3d_time_count(av) + 1 > ntdim_3d(av) ) THEN
214	WRITE ( message_string, * ) 'Output of 3d data is not given at t=', &
215	time_since_reference_point, 's because the maximum ', &
216	'number of output time levels is ', &
217	'exceeded.'
218	CALL message( 'data_output_3d', 'PA0387', 0, 1, 0, 6, 0 )
219	CALL cpu_log( log_point(14), 'data_output_3d', 'stop' )
220	RETURN
221	ENDIF
222	ENDIF
223
224	!
225	!-- Update the netCDF time axis.
226	!-- In case of parallel output, this is only done by PE0 to increase the performance.
227	#if defined( __netcdf )
228	do3d_time_count(av) = do3d_time_count(av) + 1
229	IF ( myid == 0 ) THEN
230	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_time_3d(av), &
231	(/ time_since_reference_point /), &
232	start = (/ do3d_time_count(av) /), &
233	count = (/ 1 /) )
234	CALL netcdf_handle_error( 'data_output_3d', 376 )
235	ENDIF
236	#endif
237
238	!
239	!-- Loop over all variables to be written.
240	ivar = 1
241
242	DO WHILE ( do3d(av,ivar)(1:1) /= ' ' )
243
244	!
245	!-- Initiate found flag and resorting flag
246	found = .FALSE.
247	resorted = .FALSE.
248	trimvar = TRIM( do3d(av,ivar) )
249
250	!
251	!-- Temporary solution to account for data output within the new urban surface model
252	!-- (urban_surface_mod.f90), see also SELECT CASE ( trimvar ).
253	!-- Store the array chosen on the temporary array.
254	nzb_do = nzb
255	!
256	!-- Set top index for 3D output. Note in case of plant-canopy model these index is determined by
257	!-- pch_index.
258	IF ( plant_canopy .AND. trimvar(1:4) == 'pcm_' ) THEN
259	nzt_do = pch_index
260	ELSE
261	nzt_do = nz_do3d
262	ENDIF
263
264	!
265	!-- Allocate a temporary array with the desired output dimensions.
266	ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )
267	!
268	!-- Before each output, set array local_pf to fill value
269	local_pf = REAL( fill_value, KIND=wp )
270	!
271	!-- Set masking flag for topography for not resorted arrays
272	flag_nr = 0
273
274	SELECT CASE ( trimvar )
275
276	CASE ( 'e' )
277	IF ( av == 0 ) THEN
278	to_be_resorted => e
279	ELSE
280	IF ( .NOT. ALLOCATED( e_av ) ) THEN
281	ALLOCATE( e_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
282	e_av = REAL( fill_value, KIND = wp )
283	ENDIF
284	to_be_resorted => e_av
285	ENDIF
286
287	CASE ( 'thetal' )
288	IF ( av == 0 ) THEN
289	to_be_resorted => pt
290	ELSE
291	IF ( .NOT. ALLOCATED( lpt_av ) ) THEN
292	ALLOCATE( lpt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
293	lpt_av = REAL( fill_value, KIND = wp )
294	ENDIF
295	to_be_resorted => lpt_av
296	ENDIF
297
298	CASE ( 'p' )
299	IF ( av == 0 ) THEN
300	IF ( psolver /= 'sor' ) CALL exchange_horiz( p, nbgp )
301	to_be_resorted => p
302	ELSE
303	IF ( .NOT. ALLOCATED( p_av ) ) THEN
304	ALLOCATE( p_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
305	p_av = REAL( fill_value, KIND = wp )
306	ENDIF
307	IF ( psolver /= 'sor' ) CALL exchange_horiz( p_av, nbgp )
308	to_be_resorted => p_av
309	ENDIF
310
311	CASE ( 'pc' ) ! particle concentration (requires ghostpoint exchange)
312	IF ( av == 0 ) THEN
313	IF ( time_since_reference_point >= particle_advection_start ) THEN
314	tend = prt_count
315	ELSE
316	tend = 0.0_wp
317	ENDIF
318	DO i = nxl, nxr
319	DO j = nys, nyn
320	DO k = nzb_do, nzt_do
321	local_pf(i,j,k) = tend(k,j,i)
322	ENDDO
323	ENDDO
324	ENDDO
325	resorted = .TRUE.
326	ELSE
327	IF ( .NOT. ALLOCATED( pc_av ) ) THEN
328	ALLOCATE( pc_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
329	pc_av = REAL( fill_value, KIND = wp )
330	ENDIF
331	to_be_resorted => pc_av
332	ENDIF
333
334	CASE ( 'pr' ) ! mean particle radius (effective radius)
335	IF ( av == 0 ) THEN
336	IF ( time_since_reference_point >= particle_advection_start ) THEN
337	DO i = nxl, nxr
338	DO j = nys, nyn
339	DO k = nzb_do, nzt_do
340	number_of_particles = prt_count(k,j,i)
341	IF ( number_of_particles <= 0 ) CYCLE
342	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
343	s_r2 = 0.0_wp
344	s_r3 = 0.0_wp
345	DO n = 1, number_of_particles
346	IF ( particles(n)%particle_mask ) THEN
347	s_r2 = s_r2 + particles(n)%radius*2 particles(n)%weight_factor
348	s_r3 = s_r3 + particles(n)%radius*3 particles(n)%weight_factor
349	ENDIF
350	ENDDO
351	IF ( s_r2 > 0.0_wp ) THEN
352	mean_r = s_r3 / s_r2
353	ELSE
354	mean_r = 0.0_wp
355	ENDIF
356	tend(k,j,i) = mean_r
357	ENDDO
358	ENDDO
359	ENDDO
360	ELSE
361	tend = 0.0_wp
362	ENDIF
363	DO i = nxl, nxr
364	DO j = nys, nyn
365	DO k = nzb_do, nzt_do
366	local_pf(i,j,k) = tend(k,j,i)
367	ENDDO
368	ENDDO
369	ENDDO
370	resorted = .TRUE.
371	ELSE
372	IF ( .NOT. ALLOCATED( pr_av ) ) THEN
373	ALLOCATE( pr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
374	pr_av = REAL( fill_value, KIND = wp )
375	ENDIF
376	to_be_resorted => pr_av
377	ENDIF
378
379	CASE ( 'theta' )
380	IF ( av == 0 ) THEN
381	IF ( .NOT. bulk_cloud_model ) THEN
382	to_be_resorted => pt
383	ELSE
384	DO i = nxl, nxr
385	DO j = nys, nyn
386	DO k = nzb_do, nzt_do
387	local_pf(i,j,k) = pt(k,j,i) + lv_d_cp * d_exner(k) * ql(k,j,i)
388	ENDDO
389	ENDDO
390	ENDDO
391	resorted = .TRUE.
392	ENDIF
393	ELSE
394	IF ( .NOT. ALLOCATED( pt_av ) ) THEN
395	ALLOCATE( pt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
396	pt_av = REAL( fill_value, KIND = wp )
397	ENDIF
398	to_be_resorted => pt_av
399	ENDIF
400
401	CASE ( 'q' )
402	IF ( av == 0 ) THEN
403	to_be_resorted => q
404	ELSE
405	IF ( .NOT. ALLOCATED( q_av ) ) THEN
406	ALLOCATE( q_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
407	q_av = REAL( fill_value, KIND = wp )
408	ENDIF
409	to_be_resorted => q_av
410	ENDIF
411
412	CASE ( 'ql' )
413	IF ( av == 0 ) THEN
414	to_be_resorted => ql
415	ELSE
416	IF ( .NOT. ALLOCATED( ql_av ) ) THEN
417	ALLOCATE( ql_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
418	ql_av = REAL( fill_value, KIND = wp )
419	ENDIF
420	to_be_resorted => ql_av
421	ENDIF
422
423	CASE ( 'ql_c' )
424	IF ( av == 0 ) THEN
425	to_be_resorted => ql_c
426	ELSE
427	IF ( .NOT. ALLOCATED( ql_c_av ) ) THEN
428	ALLOCATE( ql_c_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
429	ql_c_av = REAL( fill_value, KIND = wp )
430	ENDIF
431	to_be_resorted => ql_c_av
432	ENDIF
433
434	CASE ( 'ql_v' )
435	IF ( av == 0 ) THEN
436	to_be_resorted => ql_v
437	ELSE
438	IF ( .NOT. ALLOCATED( ql_v_av ) ) THEN
439	ALLOCATE( ql_v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
440	ql_v_av = REAL( fill_value, KIND = wp )
441	ENDIF
442	to_be_resorted => ql_v_av
443	ENDIF
444
445	CASE ( 'ql_vp' )
446	IF ( av == 0 ) THEN
447	IF ( time_since_reference_point >= particle_advection_start ) THEN
448	DO i = nxl, nxr
449	DO j = nys, nyn
450	DO k = nzb_do, nzt_do
451	number_of_particles = prt_count(k,j,i)
452	IF (number_of_particles <= 0) CYCLE
453	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
454	DO n = 1, number_of_particles
455	IF ( particles(n)%particle_mask ) THEN
456	tend(k,j,i) = tend(k,j,i) + &
457	particles(n)%weight_factor / prt_count(k,j,i)
458	ENDIF
459	ENDDO
460	ENDDO
461	ENDDO
462	ENDDO
463	ELSE
464	tend = 0.0_wp
465	ENDIF
466	DO i = nxl, nxr
467	DO j = nys, nyn
468	DO k = nzb_do, nzt_do
469	local_pf(i,j,k) = tend(k,j,i)
470	ENDDO
471	ENDDO
472	ENDDO
473	resorted = .TRUE.
474	ELSE
475	IF ( .NOT. ALLOCATED( ql_vp_av ) ) THEN
476	ALLOCATE( ql_vp_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
477	ql_vp_av = REAL( fill_value, KIND = wp )
478	ENDIF
479	to_be_resorted => ql_vp_av
480	ENDIF
481
482	CASE ( 'qv' )
483	IF ( av == 0 ) THEN
484	DO i = nxl, nxr
485	DO j = nys, nyn
486	DO k = nzb_do, nzt_do
487	local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
488	ENDDO
489	ENDDO
490	ENDDO
491	resorted = .TRUE.
492	ELSE
493	IF ( .NOT. ALLOCATED( qv_av ) ) THEN
494	ALLOCATE( qv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
495	qv_av = REAL( fill_value, KIND = wp )
496	ENDIF
497	to_be_resorted => qv_av
498	ENDIF
499
500	CASE ( 's' )
501	IF ( av == 0 ) THEN
502	to_be_resorted => s
503	ELSE
504	IF ( .NOT. ALLOCATED( s_av ) ) THEN
505	ALLOCATE( s_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
506	s_av = REAL( fill_value, KIND = wp )
507	ENDIF
508	to_be_resorted => s_av
509	ENDIF
510
511	CASE ( 'u' )
512	flag_nr = 1
513	IF ( av == 0 ) THEN
514	to_be_resorted => u
515	ELSE
516	IF ( .NOT. ALLOCATED( u_av ) ) THEN
517	ALLOCATE( u_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
518	u_av = REAL( fill_value, KIND = wp )
519	ENDIF
520	to_be_resorted => u_av
521	ENDIF
522
523	CASE ( 'v' )
524	flag_nr = 2
525	IF ( av == 0 ) THEN
526	to_be_resorted => v
527	ELSE
528	IF ( .NOT. ALLOCATED( v_av ) ) THEN
529	ALLOCATE( v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
530	v_av = REAL( fill_value, KIND = wp )
531	ENDIF
532	to_be_resorted => v_av
533	ENDIF
534
535	CASE ( 'thetav' )
536	IF ( av == 0 ) THEN
537	to_be_resorted => vpt
538	ELSE
539	IF ( .NOT. ALLOCATED( vpt_av ) ) THEN
540	ALLOCATE( vpt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
541	vpt_av = REAL( fill_value, KIND = wp )
542	ENDIF
543	to_be_resorted => vpt_av
544	ENDIF
545
546	CASE ( 'w' )
547	flag_nr = 3
548	IF ( av == 0 ) THEN
549	to_be_resorted => w
550	ELSE
551	IF ( .NOT. ALLOCATED( w_av ) ) THEN
552	ALLOCATE( w_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
553	w_av = REAL( fill_value, KIND = wp )
554	ENDIF
555	to_be_resorted => w_av
556	ENDIF
557
558	CASE DEFAULT
559	!
560	!-- Quantities of other modules
561	IF ( .NOT. found ) THEN
562	CALL module_interface_data_output_3d( av, trimvar, found, local_pf, fill_value, &
563	resorted, nzb_do, nzt_do &
564	)
565	ENDIF
566
567	!
568	!-- Temporary workaround: ToDo: refactor local_pf allocation
569	IF ( .NOT. found .AND. urban_surface .AND. trimvar(1:4) == 'usm_' ) THEN
570	!
571	!-- For urban model quantities, it is required to re-allocate local_pf
572	nzb_do = nz_urban_b
573	nzt_do = nz_urban_t
574
575	DEALLOCATE ( local_pf )
576	ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )
577	local_pf = fill_value
578
579	CALL usm_data_output_3d( av, trimvar, found, local_pf, nzb_do, nzt_do )
580	resorted = .TRUE.
581
582	!
583	!-- If no soil model variable was found, re-allocate local_pf
584	IF ( .NOT. found ) THEN
585	nzb_do = nzb
586	nzt_do = nz_do3d
587
588	DEALLOCATE ( local_pf )
589	ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )
590	ENDIF
591
592	ENDIF
593
594	!
595	!-- Temporary workaround: ToDo: refactor local_pf allocation
596	IF ( .NOT. found .AND. land_surface ) THEN
597	!
598	!-- For soil model quantities, it is required to re-allocate local_pf
599	nzb_do = nzb_soil
600	nzt_do = nzt_soil
601
602	DEALLOCATE ( local_pf )
603	ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )
604	local_pf = fill_value
605
606	CALL lsm_data_output_3d( av, trimvar, found, local_pf )
607	resorted = .TRUE.
608
609	!
610	!-- If no soil model variable was found, re-allocate local_pf
611	IF ( .NOT. found ) THEN
612	nzb_do = nzb
613	nzt_do = nz_do3d
614
615	DEALLOCATE ( local_pf )
616	ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do) )
617	ENDIF
618
619	ENDIF
620
621	IF ( .NOT. found ) THEN
622	message_string = 'no output available for: ' // TRIM( do3d(av,ivar) )
623	CALL message( 'data_output_3d', 'PA0182', 0, 0, 0, 6, 0 )
624	ENDIF
625
626	END SELECT
627
628	!
629	!-- Resort the array to be output, if not done above
630	IF ( .NOT. resorted ) THEN
631	DO i = nxl, nxr
632	DO j = nys, nyn
633	DO k = nzb_do, nzt_do
634	local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), &
635	REAL( fill_value, KIND = wp ), &
636	BTEST( wall_flags_total_0(k,j,i), flag_nr ) )
637	ENDDO
638	ENDDO
639	ENDDO
640	ENDIF
641
642	!
643	!-- Output of the 3D-array
644	#if defined( __parallel )
645	IF ( netcdf_data_format < 5 ) THEN
646	!
647	!-- Non-parallel netCDF output. Data is output in parallel in FORTRAN binary format here, and
648	!-- later collected into one file by combine_plot_fields
649	IF ( myid == 0 ) THEN
650	WRITE ( 30 ) time_since_reference_point, do3d_time_count(av), av
651	ENDIF
652	DO i = 0, io_blocks-1
653	IF ( i == io_group ) THEN
654	WRITE ( 30 ) nxl, nxr, nys, nyn, nzb_do, nzt_do
655	WRITE ( 30 ) REAL( local_pf(:,:,nzb_do:nzt_do), KIND = sp )
656	ENDIF
657
658	CALL MPI_BARRIER( comm2d, ierr )
659
660	ENDDO
661
662	ELSE
663	#if defined( __netcdf )
664	!
665	!-- Parallel output in netCDF4/HDF5 format.
666	! IF ( nxr == nx .AND. nyn /= ny ) THEN
667	! nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,ivar), &
668	! local_pf(nxl:nxr+1,nys:nyn,nzb_do:nzt_do), &
669	! start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /), &
670	! count = (/ nxr-nxl+2, nyn-nys+1, nzt_do-nzb_do+1, 1 /) )
671	! ELSEIF ( nxr /= nx .AND. nyn == ny ) THEN
672	! nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,ivar), &
673	! local_pf(nxl:nxr,nys:nyn+1,nzb_do:nzt_do), &
674	! start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /), &
675	! count = (/ nxr-nxl+1, nyn-nys+2, nzt_do-nzb_do+1, 1 /) )
676	! ELSEIF ( nxr == nx .AND. nyn == ny ) THEN
677	! nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,ivar), &
678	! local_pf(nxl:nxr+1,nys:nyn+1,nzb_do:nzt_do ), &
679	! start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /), &
680	! count = (/ nxr-nxl+2, nyn-nys+2, nzt_do-nzb_do+1, 1 /) )
681	! ELSE
682	!
683	!-- Call for parallel runs
684	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,ivar), &
685	local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do), &
686	start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /), &
687	count = (/ nxr-nxl+1, nyn-nys+1, nzt_do-nzb_do+1, 1 /) &
688	)
689	! ENDIF
690	CALL netcdf_handle_error( 'data_output_3d', 386 )
691	#endif
692	ENDIF
693	#else
694	#if defined( __netcdf )
695	!
696	!-- Call for non-parallel runs
697	nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,ivar), &
698	local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do), &
699	start = (/ 1, 1, 1, do3d_time_count(av) /), &
700	count = (/ nx+1, ny+1, nzt_do-nzb_do+1, 1 /) )
701	CALL netcdf_handle_error( 'data_output_3d', 446 )
702	#endif
703	#endif
704
705	ivar = ivar + 1
706
707	!
708	!-- Deallocate temporary array
709	DEALLOCATE ( local_pf )
710
711	ENDDO
712
713	CALL cpu_log( log_point(14), 'data_output_3d', 'stop' )
714
715	IF ( debug_output_timestep ) CALL debug_message( 'data_output_3d', 'end' )
716
717
718	END SUBROUTINE data_output_3d

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

Context Navigation

source: palm/trunk/SOURCE/data_output_3d.f90 @ 4867

Download in other formats: