Home

Context Navigation

source: palm/trunk/SOURCE/diagnostic_output_quantities_mod.f90 @ 4540

Last change on this file since 4540 was 4535, checked in by raasch, 4 years ago
bugfix for restart data format query
Property svn:keywords set to `Id`
File size: 79.4 KB

Line
1	!> @file diagnostic_output_quantities_mod.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
6	! terms of the GNU General Public License as published by the Free Software
7	! Foundation, either version 3 of the License, or (at your option) any later
8	! 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-2020 Leibniz Universitaet Hannover
18	!------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! ------------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: diagnostic_output_quantities_mod.f90 4535 2020-05-15 12:07:23Z raasch $
27	! bugfix for restart data format query
28	!
29	! 4518 2020-05-04 15:44:28Z suehring
30	! * Define arrays over ghost points in order to allow for standard mpi-io
31	! treatment. By this modularization of restart-data input is possible with
32	! the module interface.
33	! * Bugfix: add missing restart input of wtheta_av, wq_av, wu_av, and wv_av.
34	!
35	! 4517 2020-05-03 14:29:30Z raasch
36	! use statement for exchange horiz added,
37	! bugfix for call of exchange horiz 2d
38	!
39	! 4431 2020-02-27 23:23:01Z gronemeier
40	! added wspeed and wdir output; bugfix: set fill_value in case of masked output
41	!
42	! 4360 2020-01-07 11:25:50Z suehring
43	! added output of wu, wv, wtheta and wq to enable covariance calculation
44	! according to temporal EC method
45	!
46	! 4346 2019-12-18 11:55:56Z motisi
47	! Introduction of wall_flags_total_0, which currently sets bits based on static
48	! topography information used in wall_flags_static_0
49	!
50	! 4331 2019-12-10 18:25:02Z suehring
51	! - Modularize 2-m potential temperature output
52	! - New output for 10-m wind speed
53	!
54	! 4329 2019-12-10 15:46:36Z motisi
55	! Renamed wall_flags_0 to wall_flags_static_0
56	!
57	! 4182 2019-08-22 15:20:23Z scharf
58	! Corrected "Former revisions" section
59	!
60	! 4167 2019-08-16 11:01:48Z suehring
61	! Changed behaviour of masked output over surface to follow terrain and ignore
62	! buildings (J.Resler, T.Gronemeier)
63	!
64	! 4157 2019-08-14 09:19:12Z suehring
65	! Initialization restructured, in order to work also when data output during
66	! spin-up is enabled.
67	!
68	! 4132 2019-08-02 12:34:17Z suehring
69	! Bugfix in masked data output
70	!
71	! 4069 2019-07-01 14:05:51Z Giersch
72	! Masked output running index mid has been introduced as a local variable to
73	! avoid runtime error (Loop variable has been modified) in time_integration
74	!
75	! 4039 2019-06-18 10:32:41Z suehring
76	! - Add output of uu, vv, ww to enable variance calculation according temporal
77	! EC method
78	! - Allocate arrays only when they are required
79	! - Formatting adjustment
80	! - Rename subroutines
81	! - Further modularization
82	!
83	! 3998 2019-05-23 13:38:11Z suehring
84	! Bugfix in gathering all output strings
85	!
86	! 3995 2019-05-22 18:59:54Z suehring
87	! Avoid compiler warnings about unused variable and fix string operation which
88	! is not allowed with PGI compiler
89	!
90	! 3994 2019-05-22 18:08:09Z suehring
91	! Initial revision
92	!
93	! Authors:
94	! --------
95	! @author Farah Kanani-Suehring
96	!
97	!
98	! Description:
99	! ------------
100	!> ...
101	!------------------------------------------------------------------------------!
102	MODULE diagnostic_output_quantities_mod
103
104	USE arrays_3d, &
105	ONLY: ddzu, &
106	pt, &
107	q, &
108	u, &
109	v, &
110	w, &
111	zu, &
112	zw
113
114	USE basic_constants_and_equations_mod, &
115	ONLY: kappa, pi
116
117	USE control_parameters, &
118	ONLY: current_timestep_number, &
119	data_output, &
120	message_string, &
121	restart_data_format_output, &
122	varnamelength
123	!
124	! USE cpulog, &
125	! ONLY: cpu_log, log_point
126
127	USE exchange_horiz_mod, &
128	ONLY: exchange_horiz_2d
129
130	USE grid_variables, &
131	ONLY: ddx, ddy
132
133	USE indices, &
134	ONLY: nbgp, &
135	nxl, &
136	nxlg, &
137	nxr, &
138	nxrg, &
139	nyn, &
140	nyng, &
141	nys, &
142	nysg, &
143	nzb, &
144	nzt, &
145	wall_flags_total_0
146
147	USE kinds
148
149	USE restart_data_mpi_io_mod, &
150	ONLY: rd_mpi_io_check_array, &
151	rrd_mpi_io, &
152	wrd_mpi_io
153
154	USE surface_mod, &
155	ONLY: surf_def_h, &
156	surf_lsm_h, &
157	surf_type, &
158	surf_usm_h
159
160
161	IMPLICIT NONE
162
163	CHARACTER(LEN=varnamelength), DIMENSION(500) :: do_all = ' '
164
165	INTEGER(iwp) :: timestep_number_at_prev_calc = 0 !< ...at previous diagnostic output calculation
166
167	LOGICAL :: initialized_diagnostic_output_quantities = .FALSE. !< flag indicating whether output is initialized
168	LOGICAL :: prepared_diagnostic_output_quantities = .FALSE. !< flag indicating whether output is p
169
170	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pt_2m !< 2-m air potential temperature
171	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pt_2m_av !< averaged 2-m air potential temperature
172	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uv_10m !< horizontal wind speed at 10m
173	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uv_10m_av !< averaged horizontal wind speed at 10m
174
175	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ti !< rotation(u,v,w) aka turbulence intensity
176	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ti_av !< avg. rotation(u,v,w) aka turbulence intensity
177	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: u_center !< u at center of grid box
178	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: u_center_av !< mean of u_center
179	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: uu !< uu
180	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: uu_av !< mean of uu
181	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wspeed !< horizontal wind speed
182	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wspeed_av !< mean of horizotal wind speed
183	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: v_center !< v at center of grid box
184	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: v_center_av !< mean of v_center
185	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: vv !< vv
186	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: vv_av !< mean of vv
187	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wdir !< wind direction
188	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wdir_av !< mean wind direction
189	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ww !< ww
190	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ww_av !< mean of ww
191	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wu !< wu
192	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wu_av !< mean of wu
193	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wv !< wv
194	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wv_av !< mean of wv
195	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wtheta !< wtheta
196	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wtheta_av !< mean of wtheta
197	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wq !< wq
198	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wq_av !< mean of wq
199
200
201	SAVE
202
203	PRIVATE
204
205	!
206	!-- Public variables
207	PUBLIC do_all, &
208	initialized_diagnostic_output_quantities, &
209	prepared_diagnostic_output_quantities, &
210	timestep_number_at_prev_calc, &
211	pt_2m_av, &
212	ti_av, &
213	u_center_av, &
214	uu_av, &
215	uv_10m_av, &
216	v_center_av, &
217	vv_av, &
218	wdir_av, &
219	wspeed_av, &
220	ww_av
221	!
222	!-- Public routines
223	PUBLIC doq_3d_data_averaging, &
224	doq_calculate, &
225	doq_check_data_output, &
226	doq_define_netcdf_grid, &
227	doq_init, &
228	doq_output_2d, &
229	doq_output_3d, &
230	doq_output_mask, &
231	doq_rrd_local, &
232	doq_wrd_local
233
234
235	INTERFACE doq_3d_data_averaging
236	MODULE PROCEDURE doq_3d_data_averaging
237	END INTERFACE doq_3d_data_averaging
238
239	INTERFACE doq_calculate
240	MODULE PROCEDURE doq_calculate
241	END INTERFACE doq_calculate
242
243	INTERFACE doq_check_data_output
244	MODULE PROCEDURE doq_check_data_output
245	END INTERFACE doq_check_data_output
246
247	INTERFACE doq_define_netcdf_grid
248	MODULE PROCEDURE doq_define_netcdf_grid
249	END INTERFACE doq_define_netcdf_grid
250
251	INTERFACE doq_output_2d
252	MODULE PROCEDURE doq_output_2d
253	END INTERFACE doq_output_2d
254
255	INTERFACE doq_output_3d
256	MODULE PROCEDURE doq_output_3d
257	END INTERFACE doq_output_3d
258
259	INTERFACE doq_output_mask
260	MODULE PROCEDURE doq_output_mask
261	END INTERFACE doq_output_mask
262
263	INTERFACE doq_init
264	MODULE PROCEDURE doq_init
265	END INTERFACE doq_init
266
267	INTERFACE doq_prepare
268	MODULE PROCEDURE doq_prepare
269	END INTERFACE doq_prepare
270
271	INTERFACE doq_rrd_local
272	MODULE PROCEDURE doq_rrd_local_ftn
273	MODULE PROCEDURE doq_rrd_local_mpi
274	END INTERFACE doq_rrd_local
275
276	INTERFACE doq_wrd_local
277	MODULE PROCEDURE doq_wrd_local
278	END INTERFACE doq_wrd_local
279
280
281	CONTAINS
282
283	!------------------------------------------------------------------------------!
284	! Description:
285	! ------------
286	!> Sum up and time-average diagnostic output quantities as well as allocate
287	!> the array necessary for storing the average.
288	!------------------------------------------------------------------------------!
289	SUBROUTINE doq_3d_data_averaging( mode, variable )
290
291	USE control_parameters, &
292	ONLY: average_count_3d
293
294	CHARACTER (LEN=*) :: mode !<
295	CHARACTER (LEN=*) :: variable !<
296
297	INTEGER(iwp) :: i !<
298	INTEGER(iwp) :: j !<
299	INTEGER(iwp) :: k !<
300
301	IF ( mode == 'allocate' ) THEN
302
303	SELECT CASE ( TRIM( variable ) )
304
305	CASE ( 'ti' )
306	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
307	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
308	ENDIF
309	ti_av = 0.0_wp
310
311	CASE ( 'uu' )
312	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
313	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
314	ENDIF
315	uu_av = 0.0_wp
316
317	CASE ( 'vv' )
318	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
319	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
320	ENDIF
321	vv_av = 0.0_wp
322
323	CASE ( 'ww' )
324	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
325	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
326	ENDIF
327	ww_av = 0.0_wp
328
329	CASE ( 'wu' )
330	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
331	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
332	ENDIF
333	wu_av = 0.0_wp
334
335	CASE ( 'wv' )
336	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
337	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
338	ENDIF
339	wv_av = 0.0_wp
340
341	CASE ( 'wtheta' )
342	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
343	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
344	ENDIF
345	wtheta_av = 0.0_wp
346
347	CASE ( 'wq' )
348	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
349	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
350	ENDIF
351	wq_av = 0.0_wp
352
353	CASE ( 'theta_2m*' )
354	IF ( .NOT. ALLOCATED( pt_2m_av ) ) THEN
355	ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
356	ENDIF
357	pt_2m_av = 0.0_wp
358
359	CASE ( 'wspeed_10m*' )
360	IF ( .NOT. ALLOCATED( uv_10m_av ) ) THEN
361	ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
362	ENDIF
363	uv_10m_av = 0.0_wp
364
365	CASE ( 'wspeed' )
366	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
367	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
368	ENDIF
369	wspeed_av = 0.0_wp
370
371	CASE ( 'wdir' )
372	IF ( .NOT. ALLOCATED( u_center_av ) ) THEN
373	ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
374	ENDIF
375	IF ( .NOT. ALLOCATED( v_center_av ) ) THEN
376	ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
377	ENDIF
378	u_center_av = 0.0_wp
379	v_center_av = 0.0_wp
380
381	CASE DEFAULT
382	CONTINUE
383
384	END SELECT
385
386	ELSEIF ( mode == 'sum' ) THEN
387
388	SELECT CASE ( TRIM( variable ) )
389
390	CASE ( 'ti' )
391	IF ( ALLOCATED( ti_av ) ) THEN
392	DO i = nxl, nxr
393	DO j = nys, nyn
394	DO k = nzb, nzt+1
395	ti_av(k,j,i) = ti_av(k,j,i) + ti(k,j,i)
396	ENDDO
397	ENDDO
398	ENDDO
399	ENDIF
400
401	CASE ( 'uu' )
402	IF ( ALLOCATED( uu_av ) ) THEN
403	DO i = nxl, nxr
404	DO j = nys, nyn
405	DO k = nzb, nzt+1
406	uu_av(k,j,i) = uu_av(k,j,i) + uu(k,j,i)
407	ENDDO
408	ENDDO
409	ENDDO
410	ENDIF
411
412	CASE ( 'vv' )
413	IF ( ALLOCATED( vv_av ) ) THEN
414	DO i = nxl, nxr
415	DO j = nys, nyn
416	DO k = nzb, nzt+1
417	vv_av(k,j,i) = vv_av(k,j,i) + vv(k,j,i)
418	ENDDO
419	ENDDO
420	ENDDO
421	ENDIF
422
423	CASE ( 'ww' )
424	IF ( ALLOCATED( ww_av ) ) THEN
425	DO i = nxl, nxr
426	DO j = nys, nyn
427	DO k = nzb, nzt+1
428	ww_av(k,j,i) = ww_av(k,j,i) + ww(k,j,i)
429	ENDDO
430	ENDDO
431	ENDDO
432	ENDIF
433
434	CASE ( 'wu' )
435	IF ( ALLOCATED( wu_av ) ) THEN
436	DO i = nxl, nxr
437	DO j = nys, nyn
438	DO k = nzb, nzt+1
439	wu_av(k,j,i) = wu_av(k,j,i) + wu(k,j,i)
440	ENDDO
441	ENDDO
442	ENDDO
443	ENDIF
444
445	CASE ( 'wv' )
446	IF ( ALLOCATED( wv_av ) ) THEN
447	DO i = nxl, nxr
448	DO j = nys, nyn
449	DO k = nzb, nzt+1
450	wv_av(k,j,i) = wv_av(k,j,i) + wv(k,j,i)
451	ENDDO
452	ENDDO
453	ENDDO
454	ENDIF
455
456	CASE ( 'wtheta' )
457	IF ( ALLOCATED( wtheta_av ) ) THEN
458	DO i = nxl, nxr
459	DO j = nys, nyn
460	DO k = nzb, nzt+1
461	wtheta_av(k,j,i) = wtheta_av(k,j,i) + wtheta(k,j,i)
462	ENDDO
463	ENDDO
464	ENDDO
465	ENDIF
466
467	CASE ( 'wq' )
468	IF ( ALLOCATED( wq_av ) ) THEN
469	DO i = nxl, nxr
470	DO j = nys, nyn
471	DO k = nzb, nzt+1
472	wq_av(k,j,i) = wq_av(k,j,i) + wq(k,j,i)
473	ENDDO
474	ENDDO
475	ENDDO
476	ENDIF
477
478	CASE ( 'theta_2m*' )
479	IF ( ALLOCATED( pt_2m_av ) ) THEN
480	DO i = nxl, nxr
481	DO j = nys, nyn
482	pt_2m_av(j,i) = pt_2m_av(j,i) + pt_2m(j,i)
483	ENDDO
484	ENDDO
485	ENDIF
486
487	CASE ( 'wspeed_10m*' )
488	IF ( ALLOCATED( uv_10m_av ) ) THEN
489	DO i = nxl, nxr
490	DO j = nys, nyn
491	uv_10m_av(j,i) = uv_10m_av(j,i) + uv_10m(j,i)
492	ENDDO
493	ENDDO
494	ENDIF
495
496	CASE ( 'wspeed' )
497	IF ( ALLOCATED( wspeed_av ) ) THEN
498	DO i = nxl, nxr
499	DO j = nys, nyn
500	DO k = nzb, nzt+1
501	wspeed_av(k,j,i) = wspeed_av(k,j,i) + wspeed(k,j,i)
502	ENDDO
503	ENDDO
504	ENDDO
505	ENDIF
506
507	CASE ( 'wdir' )
508	IF ( ALLOCATED( u_center_av ) .AND. ALLOCATED( v_center_av ) ) THEN
509	DO i = nxl, nxr
510	DO j = nys, nyn
511	DO k = nzb, nzt+1
512	u_center_av(k,j,i) = u_center_av(k,j,i) + u_center(k,j,i)
513	v_center_av(k,j,i) = v_center_av(k,j,i) + v_center(k,j,i)
514	ENDDO
515	ENDDO
516	ENDDO
517	ENDIF
518
519	CASE DEFAULT
520	CONTINUE
521
522	END SELECT
523
524	ELSEIF ( mode == 'average' ) THEN
525
526	SELECT CASE ( TRIM( variable ) )
527
528	CASE ( 'ti' )
529	IF ( ALLOCATED( ti_av ) ) THEN
530	DO i = nxl, nxr
531	DO j = nys, nyn
532	DO k = nzb, nzt+1
533	ti_av(k,j,i) = ti_av(k,j,i) / REAL( average_count_3d, KIND=wp )
534	ENDDO
535	ENDDO
536	ENDDO
537	ENDIF
538
539	CASE ( 'uu' )
540	IF ( ALLOCATED( uu_av ) ) THEN
541	DO i = nxl, nxr
542	DO j = nys, nyn
543	DO k = nzb, nzt+1
544	uu_av(k,j,i) = uu_av(k,j,i) / REAL( average_count_3d, KIND=wp )
545	ENDDO
546	ENDDO
547	ENDDO
548	ENDIF
549
550	CASE ( 'vv' )
551	IF ( ALLOCATED( vv_av ) ) THEN
552	DO i = nxl, nxr
553	DO j = nys, nyn
554	DO k = nzb, nzt+1
555	vv_av(k,j,i) = vv_av(k,j,i) / REAL( average_count_3d, KIND=wp )
556	ENDDO
557	ENDDO
558	ENDDO
559	ENDIF
560
561	CASE ( 'ww' )
562	IF ( ALLOCATED( ww_av ) ) THEN
563	DO i = nxl, nxr
564	DO j = nys, nyn
565	DO k = nzb, nzt+1
566	ww_av(k,j,i) = ww_av(k,j,i) / REAL( average_count_3d, KIND=wp )
567	ENDDO
568	ENDDO
569	ENDDO
570	ENDIF
571
572	CASE ( 'wu' )
573	IF ( ALLOCATED( wu_av ) ) THEN
574	DO i = nxl, nxr
575	DO j = nys, nyn
576	DO k = nzb, nzt+1
577	wu_av(k,j,i) = wu_av(k,j,i) / REAL( average_count_3d, KIND=wp )
578	ENDDO
579	ENDDO
580	ENDDO
581	ENDIF
582
583	CASE ( 'wv' )
584	IF ( ALLOCATED( wv_av ) ) THEN
585	DO i = nxl, nxr
586	DO j = nys, nyn
587	DO k = nzb, nzt+1
588	wv_av(k,j,i) = wv_av(k,j,i) / REAL( average_count_3d, KIND=wp )
589	ENDDO
590	ENDDO
591	ENDDO
592	ENDIF
593
594	CASE ( 'wtheta' )
595	IF ( ALLOCATED( wtheta_av ) ) THEN
596	DO i = nxl, nxr
597	DO j = nys, nyn
598	DO k = nzb, nzt+1
599	wtheta_av(k,j,i) = wtheta_av(k,j,i) / REAL( average_count_3d, KIND=wp )
600	ENDDO
601	ENDDO
602	ENDDO
603	ENDIF
604
605	CASE ( 'wq' )
606	IF ( ALLOCATED( wq_av ) ) THEN
607	DO i = nxl, nxr
608	DO j = nys, nyn
609	DO k = nzb, nzt+1
610	wq_av(k,j,i) = wq_av(k,j,i) / REAL( average_count_3d, KIND=wp )
611	ENDDO
612	ENDDO
613	ENDDO
614	ENDIF
615
616	CASE ( 'theta_2m*' )
617	IF ( ALLOCATED( pt_2m_av ) ) THEN
618	DO i = nxlg, nxrg
619	DO j = nysg, nyng
620	pt_2m_av(j,i) = pt_2m_av(j,i) / REAL( average_count_3d, KIND=wp )
621	ENDDO
622	ENDDO
623	CALL exchange_horiz_2d( pt_2m_av )
624	ENDIF
625
626	CASE ( 'wspeed_10m*' )
627	IF ( ALLOCATED( uv_10m_av ) ) THEN
628	DO i = nxlg, nxrg
629	DO j = nysg, nyng
630	uv_10m_av(j,i) = uv_10m_av(j,i) / REAL( average_count_3d, KIND=wp )
631	ENDDO
632	ENDDO
633	CALL exchange_horiz_2d( uv_10m_av )
634	ENDIF
635
636	CASE ( 'wspeed' )
637	IF ( ALLOCATED( wspeed_av ) ) THEN
638	DO i = nxl, nxr
639	DO j = nys, nyn
640	DO k = nzb, nzt+1
641	wspeed_av(k,j,i) = wspeed_av(k,j,i) / REAL( average_count_3d, KIND=wp )
642	ENDDO
643	ENDDO
644	ENDDO
645	ENDIF
646
647	CASE ( 'wdir' )
648	IF ( ALLOCATED( u_center_av ) .AND. ALLOCATED( v_center_av ) ) THEN
649
650	IF ( .NOT. ALLOCATED( wdir_av ) ) THEN
651	ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
652	ENDIF
653	wdir_av = 0.0_wp
654
655	DO i = nxl, nxr
656	DO j = nys, nyn
657	DO k = nzb, nzt+1
658	u_center_av(k,j,i) = u_center_av(k,j,i) / REAL( average_count_3d, KIND=wp )
659	v_center_av(k,j,i) = v_center_av(k,j,i) / REAL( average_count_3d, KIND=wp )
660	wdir_av(k,j,i) = ATAN2( u_center_av(k,j,i), v_center_av(k,j,i) ) &
661	/ pi * 180.0_wp + 180.0_wp
662	ENDDO
663	ENDDO
664	ENDDO
665	ENDIF
666
667	END SELECT
668
669	ENDIF
670
671
672	END SUBROUTINE doq_3d_data_averaging
673
674	!------------------------------------------------------------------------------!
675	! Description:
676	! ------------
677	!> Check data output for diagnostic output
678	!------------------------------------------------------------------------------!
679	SUBROUTINE doq_check_data_output( var, unit, i, ilen, k )
680
681	IMPLICIT NONE
682
683	CHARACTER (LEN=*) :: unit !<
684	CHARACTER (LEN=*) :: var !<
685
686	INTEGER(iwp), OPTIONAL, INTENT(IN) :: i !< Current element of data_output
687	INTEGER(iwp), OPTIONAL, INTENT(IN) :: ilen !< Length of current entry in data_output
688	INTEGER(iwp), OPTIONAL, INTENT(IN) :: k !< Output is xy mode? 0 = no, 1 = yes
689
690	SELECT CASE ( TRIM( var ) )
691
692	CASE ( 'ti' )
693	unit = '1/s'
694
695	CASE ( 'uu' )
696	unit = 'm2/s2'
697
698	CASE ( 'vv' )
699	unit = 'm2/s2'
700
701	CASE ( 'ww' )
702	unit = 'm2/s2'
703
704	CASE ( 'wu' )
705	unit = 'm2/s2'
706
707	CASE ( 'wv' )
708	unit = 'm2/s2'
709
710	CASE ( 'wtheta' )
711	unit = 'Km/s'
712
713	CASE ( 'wq' )
714	unit = 'm/s'
715
716	CASE ( 'wspeed' )
717	unit = 'm/s'
718
719	CASE ( 'wdir' )
720	unit = 'degree'
721	!
722	!-- Treat horizotal cross-section output quanatities
723	CASE ( 'theta_2m', 'wspeed_10m' )
724	!
725	!-- Check if output quantity is _xy only.
726	IF ( k == 0 .OR. data_output(i)(ilen-2:ilen) /= '_xy' ) THEN
727	message_string = 'illegal value for data_output: "' // &
728	TRIM( var ) // '" & only 2d-horizontal ' // &
729	'cross sections are allowed for this value'
730	CALL message( 'diagnostic_output', 'PA0111', 1, 2, 0, 6, 0 )
731	ENDIF
732
733	IF ( TRIM( var ) == 'theta_2m*' ) unit = 'K'
734	IF ( TRIM( var ) == 'wspeed_10m*' ) unit = 'm/s'
735
736	CASE DEFAULT
737	unit = 'illegal'
738
739	END SELECT
740
741
742	END SUBROUTINE doq_check_data_output
743
744	!------------------------------------------------------------------------------!
745	!
746	! Description:
747	! ------------
748	!> Subroutine defining appropriate grid for netcdf variables.
749	!------------------------------------------------------------------------------!
750	SUBROUTINE doq_define_netcdf_grid( variable, found, grid_x, grid_y, grid_z )
751
752	IMPLICIT NONE
753
754	CHARACTER (LEN=*), INTENT(IN) :: variable !<
755	LOGICAL, INTENT(OUT) :: found !<
756	CHARACTER (LEN=*), INTENT(OUT) :: grid_x !<
757	CHARACTER (LEN=*), INTENT(OUT) :: grid_y !<
758	CHARACTER (LEN=*), INTENT(OUT) :: grid_z !<
759
760	found = .TRUE.
761
762	SELECT CASE ( TRIM( variable ) )
763	!
764	!-- s grid
765	CASE ( 'ti', 'ti_xy', 'ti_xz', 'ti_yz', &
766	'wspeed', 'wspeed_xy', 'wspeed_xz', 'wspeed_yz', &
767	'wdir', 'wdir_xy', 'wdir_xz', 'wdir_yz', &
768	'wu', 'wu_xy', 'wu_xz', 'wu_yz', &
769	'wv', 'wv_xy', 'wv_xz', 'wv_yz', &
770	'wtheta', 'wtheta_xy', 'wtheta_xz', 'wtheta_yz', &
771	'wq', 'wq_xy', 'wq_xz', 'wq_yz')
772
773	grid_x = 'x'
774	grid_y = 'y'
775	grid_z = 'zu'
776	!
777	!-- s grid surface variables
778	CASE ( 'theta_2m_xy', 'wspeed_10m' )
779
780	grid_x = 'x'
781	grid_y = 'y'
782	grid_z = 'zu'
783	!
784	!-- u grid
785	CASE ( 'uu', 'uu_xy', 'uu_xz', 'uu_yz' )
786
787	grid_x = 'xu'
788	grid_y = 'y'
789	grid_z = 'zu'
790	!
791	!-- v grid
792	CASE ( 'vv', 'vv_xy', 'vv_xz', 'vv_yz' )
793
794	grid_x = 'x'
795	grid_y = 'yv'
796	grid_z = 'zu'
797	!
798	!-- w grid
799	CASE ( 'ww', 'ww_xy', 'ww_xz', 'ww_yz' )
800
801	grid_x = 'x'
802	grid_y = 'y'
803	grid_z = 'zw'
804
805	CASE DEFAULT
806	found = .FALSE.
807	grid_x = 'none'
808	grid_y = 'none'
809	grid_z = 'none'
810
811	END SELECT
812
813
814	END SUBROUTINE doq_define_netcdf_grid
815
816	!------------------------------------------------------------------------------!
817	!
818	! Description:
819	! ------------
820	!> Subroutine defining 2D output variables
821	!------------------------------------------------------------------------------!
822	SUBROUTINE doq_output_2d( av, variable, found, grid, &
823	mode, local_pf, two_d, nzb_do, nzt_do, fill_value )
824
825
826	IMPLICIT NONE
827
828	CHARACTER (LEN=*) :: grid !<
829	CHARACTER (LEN=*) :: mode !<
830	CHARACTER (LEN=*) :: variable !<
831
832	INTEGER(iwp) :: av !< value indicating averaged or non-averaged output
833	INTEGER(iwp) :: flag_nr !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
834	INTEGER(iwp) :: i !< grid index x-direction
835	INTEGER(iwp) :: j !< grid index y-direction
836	INTEGER(iwp) :: k !< grid index z-direction
837	INTEGER(iwp) :: nzb_do !<
838	INTEGER(iwp) :: nzt_do !<
839
840	LOGICAL :: found !< true if variable is in list
841	LOGICAL :: resorted !< true if array is resorted
842	LOGICAL :: two_d !< flag parameter that indicates 2D variables (horizontal cross sections)
843
844	REAL(wp) :: fill_value !< value for the _FillValue attribute
845
846	REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !<
847	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to array which needs to be resorted for output
848
849	flag_nr = 0
850	found = .TRUE.
851	resorted = .FALSE.
852	two_d = .FALSE.
853
854	SELECT CASE ( TRIM( variable ) )
855
856	CASE ( 'ti_xy', 'ti_xz', 'ti_yz' )
857	IF ( av == 0 ) THEN
858	to_be_resorted => ti
859	ELSE
860	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
861	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
862	ti_av = REAL( fill_value, KIND = wp )
863	ENDIF
864	to_be_resorted => ti_av
865	ENDIF
866	flag_nr = 0
867
868	IF ( mode == 'xy' ) grid = 'zu'
869
870	CASE ( 'uu_xy', 'uu_xz', 'uu_yz' )
871	IF ( av == 0 ) THEN
872	to_be_resorted => uu
873	ELSE
874	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
875	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
876	uu_av = REAL( fill_value, KIND = wp )
877	ENDIF
878	to_be_resorted => uu_av
879	ENDIF
880	flag_nr = 1
881
882	IF ( mode == 'xy' ) grid = 'zu'
883
884	CASE ( 'vv_xy', 'vv_xz', 'vv_yz' )
885	IF ( av == 0 ) THEN
886	to_be_resorted => vv
887	ELSE
888	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
889	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
890	vv_av = REAL( fill_value, KIND = wp )
891	ENDIF
892	to_be_resorted => vv_av
893	ENDIF
894	flag_nr = 2
895
896	IF ( mode == 'xy' ) grid = 'zu'
897
898	CASE ( 'ww_xy', 'ww_xz', 'ww_yz' )
899	IF ( av == 0 ) THEN
900	to_be_resorted => ww
901	ELSE
902	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
903	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
904	ww_av = REAL( fill_value, KIND = wp )
905	ENDIF
906	to_be_resorted => ww_av
907	ENDIF
908	flag_nr = 3
909
910	IF ( mode == 'xy' ) grid = 'zw'
911
912	CASE ( 'wu_xy', 'wu_xz', 'wu_yz' )
913	IF ( av == 0 ) THEN
914	to_be_resorted => wu
915	ELSE
916	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
917	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
918	wu_av = REAL( fill_value, KIND = wp )
919	ENDIF
920	to_be_resorted => wu_av
921	ENDIF
922	flag_nr = 0
923
924	IF ( mode == 'xy' ) grid = 'zw'
925
926	CASE ( 'wv_xy', 'wv_xz', 'wv_yz' )
927	IF ( av == 0 ) THEN
928	to_be_resorted => wv
929	ELSE
930	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
931	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
932	wv_av = REAL( fill_value, KIND = wp )
933	ENDIF
934	to_be_resorted => wv_av
935	ENDIF
936	flag_nr = 0
937
938	IF ( mode == 'xy' ) grid = 'zw'
939
940	CASE ( 'wtheta_xy', 'wtheta_xz', 'wtheta_yz' )
941	IF ( av == 0 ) THEN
942	to_be_resorted => wtheta
943	ELSE
944	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
945	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
946	wtheta_av = REAL( fill_value, KIND = wp )
947	ENDIF
948	to_be_resorted => wtheta_av
949	ENDIF
950	flag_nr = 0
951
952	IF ( mode == 'xy' ) grid = 'zw'
953
954	CASE ( 'wq_xy', 'wq_xz', 'wq_yz' )
955	IF ( av == 0 ) THEN
956	to_be_resorted => wq
957	ELSE
958	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
959	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
960	wq_av = REAL( fill_value, KIND = wp )
961	ENDIF
962	to_be_resorted => wq_av
963	ENDIF
964	flag_nr = 0
965
966	IF ( mode == 'xy' ) grid = 'zw'
967
968	CASE ( 'theta_2m*_xy' ) ! 2d-array
969	IF ( av == 0 ) THEN
970	DO i = nxl, nxr
971	DO j = nys, nyn
972	local_pf(i,j,nzb+1) = pt_2m(j,i)
973	ENDDO
974	ENDDO
975	ELSE
976	IF ( .NOT. ALLOCATED( pt_2m_av ) ) THEN
977	ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
978	pt_2m_av = REAL( fill_value, KIND = wp )
979	ENDIF
980	DO i = nxl, nxr
981	DO j = nys, nyn
982	local_pf(i,j,nzb+1) = pt_2m_av(j,i)
983	ENDDO
984	ENDDO
985	ENDIF
986	resorted = .TRUE.
987	two_d = .TRUE.
988	grid = 'zu1'
989
990	CASE ( 'wspeed_10m*_xy' ) ! 2d-array
991	IF ( av == 0 ) THEN
992	DO i = nxl, nxr
993	DO j = nys, nyn
994	local_pf(i,j,nzb+1) = uv_10m(j,i)
995	ENDDO
996	ENDDO
997	ELSE
998	IF ( .NOT. ALLOCATED( uv_10m_av ) ) THEN
999	ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
1000	uv_10m_av = REAL( fill_value, KIND = wp )
1001	ENDIF
1002	DO i = nxl, nxr
1003	DO j = nys, nyn
1004	local_pf(i,j,nzb+1) = uv_10m_av(j,i)
1005	ENDDO
1006	ENDDO
1007	ENDIF
1008	resorted = .TRUE.
1009	two_d = .TRUE.
1010	grid = 'zu1'
1011
1012	CASE ( 'wspeed_xy', 'wspeed_xz', 'wspeed_yz' )
1013	IF ( av == 0 ) THEN
1014	to_be_resorted => wspeed
1015	ELSE
1016	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
1017	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1018	wspeed_av = REAL( fill_value, KIND = wp )
1019	ENDIF
1020	to_be_resorted => wspeed_av
1021	ENDIF
1022	flag_nr = 0
1023
1024	IF ( mode == 'xy' ) grid = 'zu'
1025
1026	CASE ( 'wdir_xy', 'wdir_xz', 'wdir_yz' )
1027	IF ( av == 0 ) THEN
1028	to_be_resorted => wdir
1029	ELSE
1030	IF ( .NOT. ALLOCATED( wdir_av ) ) THEN
1031	ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1032	wdir_av = REAL( fill_value, KIND = wp )
1033	ENDIF
1034	to_be_resorted => wdir_av
1035	ENDIF
1036	flag_nr = 0
1037
1038	IF ( mode == 'xy' ) grid = 'zu'
1039
1040	CASE DEFAULT
1041	found = .FALSE.
1042	grid = 'none'
1043
1044	END SELECT
1045
1046	IF ( found .AND. .NOT. resorted ) THEN
1047	DO i = nxl, nxr
1048	DO j = nys, nyn
1049	DO k = nzb_do, nzt_do
1050	local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), &
1051	REAL( fill_value, KIND = wp ), &
1052	BTEST( wall_flags_total_0(k,j,i), flag_nr ) )
1053	ENDDO
1054	ENDDO
1055	ENDDO
1056	ENDIF
1057
1058	END SUBROUTINE doq_output_2d
1059
1060
1061	!------------------------------------------------------------------------------!
1062	!
1063	! Description:
1064	! ------------
1065	!> Subroutine defining 3D output variables
1066	!------------------------------------------------------------------------------!
1067	SUBROUTINE doq_output_3d( av, variable, found, local_pf, fill_value, nzb_do, &
1068	nzt_do )
1069
1070	IMPLICIT NONE
1071
1072	CHARACTER (LEN=*) :: variable !<
1073
1074	INTEGER(iwp) :: av !< index indicating averaged or instantaneous output
1075	INTEGER(iwp) :: flag_nr !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
1076	INTEGER(iwp) :: i !< index variable along x-direction
1077	INTEGER(iwp) :: j !< index variable along y-direction
1078	INTEGER(iwp) :: k !< index variable along z-direction
1079	INTEGER(iwp) :: nzb_do !< lower limit of the data output (usually 0)
1080	INTEGER(iwp) :: nzt_do !< vertical upper limit of the data output (usually nz_do3d)
1081
1082	LOGICAL :: found !< true if variable is in list
1083	LOGICAL :: resorted !< true if array is resorted
1084
1085	REAL(wp) :: fill_value !< value for the _FillValue attribute
1086
1087	REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !<
1088	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to array which needs to be resorted for output
1089
1090	flag_nr = 0
1091	found = .TRUE.
1092	resorted = .FALSE.
1093
1094	SELECT CASE ( TRIM( variable ) )
1095
1096	CASE ( 'ti' )
1097	IF ( av == 0 ) THEN
1098	to_be_resorted => ti
1099	ELSE
1100	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
1101	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1102	ti_av = REAL( fill_value, KIND = wp )
1103	ENDIF
1104	to_be_resorted => ti_av
1105	ENDIF
1106	flag_nr = 0
1107
1108	CASE ( 'uu' )
1109	IF ( av == 0 ) THEN
1110	to_be_resorted => uu
1111	ELSE
1112	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
1113	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1114	uu_av = REAL( fill_value, KIND = wp )
1115	ENDIF
1116	to_be_resorted => uu_av
1117	ENDIF
1118	flag_nr = 1
1119
1120	CASE ( 'vv' )
1121	IF ( av == 0 ) THEN
1122	to_be_resorted => vv
1123	ELSE
1124	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
1125	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1126	vv_av = REAL( fill_value, KIND = wp )
1127	ENDIF
1128	to_be_resorted => vv_av
1129	ENDIF
1130	flag_nr = 2
1131
1132	CASE ( 'ww' )
1133	IF ( av == 0 ) THEN
1134	to_be_resorted => ww
1135	ELSE
1136	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
1137	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1138	ww_av = REAL( fill_value, KIND = wp )
1139	ENDIF
1140	to_be_resorted => ww_av
1141	ENDIF
1142	flag_nr = 3
1143
1144	CASE ( 'wu' )
1145	IF ( av == 0 ) THEN
1146	to_be_resorted => wu
1147	ELSE
1148	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
1149	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1150	wu_av = REAL( fill_value, KIND = wp )
1151	ENDIF
1152	to_be_resorted => wu_av
1153	ENDIF
1154	flag_nr = 0
1155
1156	CASE ( 'wv' )
1157	IF ( av == 0 ) THEN
1158	to_be_resorted => wv
1159	ELSE
1160	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
1161	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1162	wv_av = REAL( fill_value, KIND = wp )
1163	ENDIF
1164	to_be_resorted => wv_av
1165	ENDIF
1166	flag_nr = 0
1167
1168	CASE ( 'wtheta' )
1169	IF ( av == 0 ) THEN
1170	to_be_resorted => wtheta
1171	ELSE
1172	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
1173	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1174	wtheta_av = REAL( fill_value, KIND = wp )
1175	ENDIF
1176	to_be_resorted => wtheta_av
1177	ENDIF
1178	flag_nr = 0
1179
1180	CASE ( 'wq' )
1181	IF ( av == 0 ) THEN
1182	to_be_resorted => wq
1183	ELSE
1184	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
1185	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1186	wq_av = REAL( fill_value, KIND = wp )
1187	ENDIF
1188	to_be_resorted => wq_av
1189	ENDIF
1190	flag_nr = 0
1191
1192	CASE ( 'wspeed' )
1193	IF ( av == 0 ) THEN
1194	to_be_resorted => wspeed
1195	ELSE
1196	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
1197	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1198	wspeed_av = REAL( fill_value, KIND = wp )
1199	ENDIF
1200	to_be_resorted => wspeed_av
1201	ENDIF
1202	flag_nr = 0
1203
1204	CASE ( 'wdir' )
1205	IF ( av == 0 ) THEN
1206	to_be_resorted => wdir
1207	ELSE
1208	IF ( .NOT. ALLOCATED( wdir_av ) ) THEN
1209	ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1210	wdir_av = REAL( fill_value, KIND = wp )
1211	ENDIF
1212	to_be_resorted => wdir_av
1213	ENDIF
1214	flag_nr = 0
1215
1216	CASE DEFAULT
1217	found = .FALSE.
1218
1219	END SELECT
1220
1221	IF ( found .AND. .NOT. resorted ) THEN
1222	DO i = nxl, nxr
1223	DO j = nys, nyn
1224	DO k = nzb_do, nzt_do
1225	local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), &
1226	REAL( fill_value, KIND = wp ), &
1227	BTEST( wall_flags_total_0(k,j,i), flag_nr ) )
1228	ENDDO
1229	ENDDO
1230	ENDDO
1231	ENDIF
1232
1233	END SUBROUTINE doq_output_3d
1234
1235	! Description:
1236	! ------------
1237	!> Resorts the user-defined output quantity with indices (k,j,i) to a
1238	!> temporary array with indices (i,j,k) for masked data output.
1239	!------------------------------------------------------------------------------!
1240	SUBROUTINE doq_output_mask( av, variable, found, local_pf, mid )
1241
1242	USE control_parameters
1243
1244	USE indices
1245
1246	IMPLICIT NONE
1247
1248	CHARACTER (LEN=*) :: variable !<
1249	CHARACTER (LEN=5) :: grid !< flag to distinquish between staggered grids
1250
1251	INTEGER(iwp) :: av !< index indicating averaged or instantaneous output
1252	INTEGER(iwp) :: flag_nr !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
1253	INTEGER(iwp) :: i !< index variable along x-direction
1254	INTEGER(iwp) :: j !< index variable along y-direction
1255	INTEGER(iwp) :: k !< index variable along z-direction
1256	INTEGER(iwp) :: im !< loop index for masked variables
1257	INTEGER(iwp) :: jm !< loop index for masked variables
1258	INTEGER(iwp) :: kk !< masked output index variable along z-direction
1259	INTEGER(iwp) :: mid !< masked output running index
1260	INTEGER(iwp) :: ktt !< k index of highest horizontal surface
1261
1262	LOGICAL :: found !< true if variable is in list
1263	LOGICAL :: resorted !< true if array is resorted
1264
1265	REAL(wp), &
1266	DIMENSION(mask_size_l(mid,1),mask_size_l(mid,2),mask_size_l(mid,3)) :: &
1267	local_pf !<
1268	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to array which needs to be resorted for output
1269
1270	REAL(wp), PARAMETER :: fill_value = -9999.0_wp !< value for the _FillValue attribute
1271
1272	flag_nr = 0
1273	found = .TRUE.
1274	resorted = .FALSE.
1275	grid = 's'
1276
1277	SELECT CASE ( TRIM( variable ) )
1278
1279	CASE ( 'ti' )
1280	IF ( av == 0 ) THEN
1281	to_be_resorted => ti
1282	ELSE
1283	to_be_resorted => ti_av
1284	ENDIF
1285	grid = 's'
1286	flag_nr = 0
1287
1288	CASE ( 'uu' )
1289	IF ( av == 0 ) THEN
1290	to_be_resorted => uu
1291	ELSE
1292	to_be_resorted => uu_av
1293	ENDIF
1294	grid = 'u'
1295	flag_nr = 1
1296
1297	CASE ( 'vv' )
1298	IF ( av == 0 ) THEN
1299	to_be_resorted => vv
1300	ELSE
1301	to_be_resorted => vv_av
1302	ENDIF
1303	grid = 'v'
1304	flag_nr = 2
1305
1306	CASE ( 'ww' )
1307	IF ( av == 0 ) THEN
1308	to_be_resorted => ww
1309	ELSE
1310	to_be_resorted => ww_av
1311	ENDIF
1312	grid = 'w'
1313	flag_nr = 3
1314
1315	CASE ( 'wu' )
1316	IF ( av == 0 ) THEN
1317	to_be_resorted => wu
1318	ELSE
1319	to_be_resorted => wu_av
1320	ENDIF
1321	grid = 's'
1322	flag_nr = 0
1323
1324	CASE ( 'wv' )
1325	IF ( av == 0 ) THEN
1326	to_be_resorted => wv
1327	ELSE
1328	to_be_resorted => wv_av
1329	ENDIF
1330	grid = 's'
1331	flag_nr = 0
1332
1333	CASE ( 'wtheta' )
1334	IF ( av == 0 ) THEN
1335	to_be_resorted => wtheta
1336	ELSE
1337	to_be_resorted => wtheta_av
1338	ENDIF
1339	grid = 's'
1340	flag_nr = 0
1341
1342	CASE ( 'wq' )
1343	IF ( av == 0 ) THEN
1344	to_be_resorted => wq
1345	ELSE
1346	to_be_resorted => wq_av
1347	ENDIF
1348	grid = 's'
1349	flag_nr = 0
1350
1351	CASE ( 'wspeed' )
1352	IF ( av == 0 ) THEN
1353	to_be_resorted => wspeed
1354	ELSE
1355	to_be_resorted => wspeed_av
1356	ENDIF
1357	grid = 's'
1358	flag_nr = 0
1359
1360	CASE ( 'wdir' )
1361	IF ( av == 0 ) THEN
1362	to_be_resorted => wdir
1363	ELSE
1364	to_be_resorted => wdir_av
1365	ENDIF
1366	grid = 's'
1367	flag_nr = 0
1368
1369	CASE DEFAULT
1370	found = .FALSE.
1371
1372	END SELECT
1373
1374	IF ( found .AND. .NOT. resorted ) THEN
1375	IF ( .NOT. mask_surface(mid) ) THEN
1376	!
1377	!-- Default masked output
1378	DO i = 1, mask_size_l(mid,1)
1379	DO j = 1, mask_size_l(mid,2)
1380	DO k = 1, mask_size_l(mid,3)
1381	local_pf(i,j,k) = MERGE( to_be_resorted(mask_k(mid,k), &
1382	mask_j(mid,j), &
1383	mask_i(mid,i)), &
1384	REAL( fill_value, KIND = wp ), &
1385	BTEST( wall_flags_total_0( &
1386	mask_k(mid,k), &
1387	mask_j(mid,j), &
1388	mask_i(mid,i)), &
1389	flag_nr ) )
1390	ENDDO
1391	ENDDO
1392	ENDDO
1393
1394	ELSE
1395	!
1396	!-- Terrain-following masked output
1397	DO i = 1, mask_size_l(mid,1)
1398	DO j = 1, mask_size_l(mid,2)
1399	!
1400	!-- Get k index of the highest terraing surface
1401	im = mask_i(mid,i)
1402	jm = mask_j(mid,j)
1403	ktt = MINLOC( MERGE( 1, 0, BTEST( wall_flags_total_0(:,jm,im), 5 )), &
1404	DIM = 1 ) - 1
1405	DO k = 1, mask_size_l(mid,3)
1406	kk = MIN( ktt+mask_k(mid,k), nzt+1 )
1407	!
1408	!-- Set value if not in building
1409	IF ( BTEST( wall_flags_total_0(kk,jm,im), 6 ) ) THEN
1410	local_pf(i,j,k) = fill_value
1411	ELSE
1412	local_pf(i,j,k) = to_be_resorted(kk,jm,im)
1413	ENDIF
1414	ENDDO
1415	ENDDO
1416	ENDDO
1417
1418	ENDIF
1419	ENDIF
1420
1421	END SUBROUTINE doq_output_mask
1422
1423	!------------------------------------------------------------------------------!
1424	! Description:
1425	! ------------
1426	!> Allocate required arrays
1427	!------------------------------------------------------------------------------!
1428	SUBROUTINE doq_init
1429
1430	IMPLICIT NONE
1431
1432	INTEGER(iwp) :: ivar !< loop index over all 2d/3d/mask output quantities
1433
1434	!
1435	!-- Next line is to avoid compiler warnings about unused variables
1436	IF ( timestep_number_at_prev_calc == 0 ) CONTINUE
1437	!
1438	!-- Preparatory steps and initialization of output arrays
1439	IF ( .NOT. prepared_diagnostic_output_quantities ) CALL doq_prepare
1440
1441	initialized_diagnostic_output_quantities = .FALSE.
1442
1443	ivar = 1
1444
1445	DO WHILE ( ivar <= SIZE( do_all ) )
1446
1447	SELECT CASE ( TRIM( do_all(ivar) ) )
1448	!
1449	!-- Allocate array for 'turbulence intensity'
1450	CASE ( 'ti' )
1451	IF ( .NOT. ALLOCATED( ti ) ) THEN
1452	ALLOCATE( ti(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1453	ti = 0.0_wp
1454	ENDIF
1455	!
1456	!-- Allocate array for uu
1457	CASE ( 'uu' )
1458	IF ( .NOT. ALLOCATED( uu ) ) THEN
1459	ALLOCATE( uu(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1460	uu = 0.0_wp
1461	ENDIF
1462	!
1463	!-- Allocate array for vv
1464	CASE ( 'vv' )
1465	IF ( .NOT. ALLOCATED( vv ) ) THEN
1466	ALLOCATE( vv(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1467	vv = 0.0_wp
1468	ENDIF
1469	!
1470	!-- Allocate array for ww
1471	CASE ( 'ww' )
1472	IF ( .NOT. ALLOCATED( ww ) ) THEN
1473	ALLOCATE( ww(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1474	ww = 0.0_wp
1475	ENDIF
1476	!
1477	!-- Allocate array for wu
1478	CASE ( 'wu' )
1479	IF ( .NOT. ALLOCATED( wu ) ) THEN
1480	ALLOCATE( wu(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1481	wu = 0.0_wp
1482	ENDIF
1483	!
1484	!-- Allocate array for wv
1485	CASE ( 'wv' )
1486	IF ( .NOT. ALLOCATED( wv ) ) THEN
1487	ALLOCATE( wv(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1488	wv = 0.0_wp
1489	ENDIF
1490	!
1491	!-- Allocate array for wtheta
1492	CASE ( 'wtheta' )
1493	IF ( .NOT. ALLOCATED( wtheta ) ) THEN
1494	ALLOCATE( wtheta(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1495	wtheta = 0.0_wp
1496	ENDIF
1497	!
1498	!-- Allocate array for wq
1499	CASE ( 'wq' )
1500	IF ( .NOT. ALLOCATED( wq ) ) THEN
1501	ALLOCATE( wq(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1502	wq = 0.0_wp
1503	ENDIF
1504	!
1505	!-- Allocate array for 2-m potential temperature
1506	CASE ( 'theta_2m*' )
1507	IF ( .NOT. ALLOCATED( pt_2m ) ) THEN
1508	ALLOCATE( pt_2m(nysg:nyng,nxlg:nxrg) )
1509	pt_2m = 0.0_wp
1510	ENDIF
1511	!
1512	!-- Allocate array for 10-m wind speed
1513	CASE ( 'wspeed_10m*' )
1514	IF ( .NOT. ALLOCATED( uv_10m ) ) THEN
1515	ALLOCATE( uv_10m(nysg:nyng,nxlg:nxrg) )
1516	uv_10m = 0.0_wp
1517	ENDIF
1518	!
1519	!-- Allocate array for wspeed
1520	CASE ( 'wspeed' )
1521	IF ( .NOT. ALLOCATED( wspeed ) ) THEN
1522	ALLOCATE( wspeed(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1523	wspeed = 0.0_wp
1524	ENDIF
1525
1526	!
1527	!-- Allocate array for wdir
1528	CASE ( 'wdir' )
1529	IF ( .NOT. ALLOCATED( u_center ) ) THEN
1530	ALLOCATE( u_center(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1531	u_center = 0.0_wp
1532	ENDIF
1533	IF ( .NOT. ALLOCATED( v_center ) ) THEN
1534	ALLOCATE( v_center(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1535	v_center = 0.0_wp
1536	ENDIF
1537	IF ( .NOT. ALLOCATED( wdir ) ) THEN
1538	ALLOCATE( wdir(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1539	wdir = 0.0_wp
1540	ENDIF
1541
1542	END SELECT
1543
1544	ivar = ivar + 1
1545	ENDDO
1546
1547	initialized_diagnostic_output_quantities = .TRUE.
1548
1549	END SUBROUTINE doq_init
1550
1551
1552	!--------------------------------------------------------------------------------------------------!
1553	! Description:
1554	! ------------
1555	!> Calculate diagnostic quantities
1556	!--------------------------------------------------------------------------------------------------!
1557	SUBROUTINE doq_calculate
1558
1559	IMPLICIT NONE
1560
1561	INTEGER(iwp) :: i !< grid index x-dimension
1562	INTEGER(iwp) :: j !< grid index y-dimension
1563	INTEGER(iwp) :: k !< grid index z-dimension
1564	INTEGER(iwp) :: ivar !< loop index over all 2d/3d/mask output quantities
1565
1566	TYPE(surf_type), POINTER :: surf !< surf-type array, used to generalize subroutines
1567
1568
1569	! CALL cpu_log( log_point(41), 'calculate_quantities', 'start' )
1570
1571	!
1572	!-- Save timestep number to check in time_integration if doq_calculate
1573	!-- has been called already, since the CALL occurs at two locations, but the calculations need to be
1574	!-- done only once per timestep.
1575	timestep_number_at_prev_calc = current_timestep_number
1576
1577	ivar = 1
1578
1579	DO WHILE ( ivar <= SIZE( do_all ) )
1580
1581	SELECT CASE ( TRIM( do_all(ivar) ) )
1582	!
1583	!-- Calculate 'turbulence intensity' from rot[(u,v,w)] at scalar grid point
1584	CASE ( 'ti' )
1585	DO i = nxl, nxr
1586	DO j = nys, nyn
1587	DO k = nzb+1, nzt
1588	ti(k,j,i) = 0.25_wp * SQRT( &
1589	( ( w(k,j+1,i) + w(k-1,j+1,i) &
1590	- w(k,j-1,i) - w(k-1,j-1,i) ) * ddy &
1591	- ( v(k+1,j,i) + v(k+1,j+1,i) &
1592	- v(k-1,j,i) - v(k-1,j+1,i) ) * ddzu(k) )**2 &
1593	+ ( ( u(k+1,j,i) + u(k+1,j,i+1) &
1594	- u(k-1,j,i) - u(k-1,j,i+1) ) * ddzu(k) &
1595	- ( w(k,j,i+1) + w(k-1,j,i+1) &
1596	- w(k,j,i-1) - w(k-1,j,i-1) ) * ddx )**2 &
1597	+ ( ( v(k,j,i+1) + v(k,j+1,i+1) &
1598	- v(k,j,i-1) - v(k,j+1,i-1) ) * ddx &
1599	- ( u(k,j+1,i) + u(k,j+1,i+1) &
1600	- u(k,j-1,i) - u(k,j-1,i+1) ) * ddy )**2 ) &
1601	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0) )
1602	ENDDO
1603	ENDDO
1604	ENDDO
1605	!
1606	!-- uu
1607	CASE ( 'uu' )
1608	DO i = nxl, nxr
1609	DO j = nys, nyn
1610	DO k = nzb+1, nzt
1611	uu(k,j,i) = u(k,j,i) * u(k,j,i) &
1612	* MERGE( 1.0_wp, 0.0_wp, &
1613	BTEST( wall_flags_total_0(k,j,i), 1) )
1614	ENDDO
1615	ENDDO
1616	ENDDO
1617	!
1618	!-- vv
1619	CASE ( 'vv' )
1620	DO i = nxl, nxr
1621	DO j = nys, nyn
1622	DO k = nzb+1, nzt
1623	vv(k,j,i) = v(k,j,i) * v(k,j,i) &
1624	* MERGE( 1.0_wp, 0.0_wp, &
1625	BTEST( wall_flags_total_0(k,j,i), 2) )
1626	ENDDO
1627	ENDDO
1628	ENDDO
1629	!
1630	!-- ww
1631	CASE ( 'ww' )
1632	DO i = nxl, nxr
1633	DO j = nys, nyn
1634	DO k = nzb+1, nzt-1
1635	ww(k,j,i) = w(k,j,i) * w(k,j,i) &
1636	* MERGE( 1.0_wp, 0.0_wp, &
1637	BTEST( wall_flags_total_0(k,j,i), 3) )
1638	ENDDO
1639	ENDDO
1640	ENDDO
1641	!
1642	!-- wu
1643	CASE ( 'wu' )
1644	DO i = nxl, nxr
1645	DO j = nys, nyn
1646	DO k = nzb+1, nzt-1
1647	wu(k,j,i) = w(k,j,i) * u(k,j,i) &
1648	* MERGE( 1.0_wp, 0.0_wp, &
1649	BTEST( wall_flags_total_0(k,j,i), 0) )
1650	ENDDO
1651	ENDDO
1652	ENDDO
1653	!
1654	!-- wv
1655	CASE ( 'wv' )
1656	DO i = nxl, nxr
1657	DO j = nys, nyn
1658	DO k = nzb+1, nzt-1
1659	wv(k,j,i) = w(k,j,i) * v(k,j,i) &
1660	* MERGE( 1.0_wp, 0.0_wp, &
1661	BTEST( wall_flags_total_0(k,j,i), 0) )
1662	ENDDO
1663	ENDDO
1664	ENDDO
1665	!
1666	!-- wtheta
1667	CASE ( 'wtheta' )
1668	DO i = nxl, nxr
1669	DO j = nys, nyn
1670	DO k = nzb+1, nzt-1
1671	wtheta(k,j,i) = w(k,j,i) * pt(k,j,i) &
1672	* MERGE( 1.0_wp, 0.0_wp, &
1673	BTEST( wall_flags_total_0(k,j,i), 0) )
1674	ENDDO
1675	ENDDO
1676	ENDDO
1677	!
1678	!-- wq
1679	CASE ( 'wq' )
1680	DO i = nxl, nxr
1681	DO j = nys, nyn
1682	DO k = nzb+1, nzt-1
1683	wq(k,j,i) = w(k,j,i) * q(k,j,i) &
1684	* MERGE( 1.0_wp, 0.0_wp, &
1685	BTEST( wall_flags_total_0(k,j,i), 0) )
1686	ENDDO
1687	ENDDO
1688	ENDDO
1689	!
1690	!-- 2-m potential temperature
1691	CASE ( 'theta_2m*' )
1692	!
1693	!-- 2-m potential temperature is caluclated from surface arrays. In
1694	!-- case the 2m level is below the first grid point, MOST is applied,
1695	!-- else, linear interpolation between two vertical grid levels is
1696	!-- applied. To access all surfaces, iterate over all horizontally-
1697	!-- upward facing surface types.
1698	surf => surf_def_h(0)
1699	CALL calc_pt_2m
1700	surf => surf_lsm_h
1701	CALL calc_pt_2m
1702	surf => surf_usm_h
1703	CALL calc_pt_2m
1704	!
1705	!-- 10-m wind speed
1706	CASE ( 'wspeed_10m*' )
1707	!
1708	!-- 10-m wind speed is caluclated from surface arrays. In
1709	!-- case the 10m level is below the first grid point, MOST is applied,
1710	!-- else, linear interpolation between two vertical grid levels is
1711	!-- applied. To access all surfaces, iterate over all horizontally-
1712	!-- upward facing surface types.
1713	surf => surf_def_h(0)
1714	CALL calc_wind_10m
1715	surf => surf_lsm_h
1716	CALL calc_wind_10m
1717	surf => surf_usm_h
1718	CALL calc_wind_10m
1719	!
1720	!-- horizontal wind speed
1721	CASE ( 'wspeed' )
1722	DO i = nxl, nxr
1723	DO j = nys, nyn
1724	DO k = nzb, nzt+1
1725	wspeed(k,j,i) = SQRT( ( 0.5_wp * ( u(k,j,i) + u(k,j,i+1) ) )**2 &
1726	+ ( 0.5_wp * ( v(k,j,i) + v(k,j+1,i) ) )**2 ) &
1727	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0) )
1728	ENDDO
1729	ENDDO
1730	ENDDO
1731
1732	!
1733	!-- horizontal wind direction
1734	CASE ( 'wdir' )
1735	DO i = nxl, nxr
1736	DO j = nys, nyn
1737	DO k = nzb, nzt+1
1738	u_center(k,j,i) = 0.5_wp * ( u(k,j,i) + u(k,j,i+1) )
1739	v_center(k,j,i) = 0.5_wp * ( v(k,j,i) + v(k,j+1,i) )
1740
1741	wdir(k,j,i) = ATAN2( u_center(k,j,i), v_center(k,j,i) ) &
1742	/ pi * 180.0_wp + 180.0_wp
1743	ENDDO
1744	ENDDO
1745	ENDDO
1746
1747	END SELECT
1748
1749	ivar = ivar + 1
1750	ENDDO
1751
1752	! CALL cpu_log( log_point(41), 'calculate_quantities', 'stop' )
1753
1754	!
1755	!-- The following block contains subroutines to calculate diagnostic
1756	!-- surface quantities.
1757	CONTAINS
1758	!------------------------------------------------------------------------------!
1759	! Description:
1760	! ------------
1761	!> Calculation of 2-m potential temperature.
1762	!------------------------------------------------------------------------------!
1763	SUBROUTINE calc_pt_2m
1764
1765	USE surface_layer_fluxes_mod, &
1766	ONLY: psi_h
1767
1768	IMPLICIT NONE
1769
1770	INTEGER(iwp) :: kk !< running index along the z-dimension
1771	INTEGER(iwp) :: m !< running index for surface elements
1772
1773	DO m = 1, surf%ns
1774
1775	i = surf%i(m)
1776	j = surf%j(m)
1777	k = surf%k(m)
1778	!
1779	!-- If 2-m level is below the first grid level, MOST is
1780	!-- used for calculation of 2-m temperature.
1781	IF ( surf%z_mo(m) > 2.0_wp ) THEN
1782	pt_2m(j,i) = surf%pt_surface(m) + surf%ts(m) / kappa &
1783	* ( LOG( 2.0_wp / surf%z0h(m) ) &
1784	- psi_h( 2.0_wp / surf%ol(m) ) &
1785	+ psi_h( surf%z0h(m) / surf%ol(m) ) )
1786	!
1787	!-- If 2-m level is above the first grid level, 2-m temperature
1788	!-- is linearly interpolated between the two nearest vertical grid
1789	!-- levels. Note, since 2-m temperature is only computed for
1790	!-- horizontal upward-facing surfaces, only a vertical
1791	!-- interpolation is necessary.
1792	ELSE
1793	!
1794	!-- zw(k-1) defines the height of the surface.
1795	kk = k
1796	DO WHILE ( zu(kk) - zw(k-1) < 2.0_wp .AND. kk <= nzt )
1797	kk = kk + 1
1798	ENDDO
1799	!
1800	!-- kk defines the index of the first grid level >= 2m.
1801	pt_2m(j,i) = pt(kk-1,j,i) + &
1802	( zw(k-1) + 2.0_wp - zu(kk-1) ) * &
1803	( pt(kk,j,i) - pt(kk-1,j,i) ) / &
1804	( zu(kk) - zu(kk-1) )
1805	ENDIF
1806
1807	ENDDO
1808
1809	END SUBROUTINE calc_pt_2m
1810
1811	!------------------------------------------------------------------------------!
1812	! Description:
1813	! ------------
1814	!> Calculation of 10-m wind speed.
1815	!------------------------------------------------------------------------------!
1816	SUBROUTINE calc_wind_10m
1817
1818	USE surface_layer_fluxes_mod, &
1819	ONLY: psi_m
1820
1821	IMPLICIT NONE
1822
1823	INTEGER(iwp) :: kk !< running index along the z-dimension
1824	INTEGER(iwp) :: m !< running index for surface elements
1825
1826	REAL(wp) :: uv_l !< wind speed at lower grid point
1827	REAL(wp) :: uv_u !< wind speed at upper grid point
1828
1829	DO m = 1, surf%ns
1830
1831	i = surf%i(m)
1832	j = surf%j(m)
1833	k = surf%k(m)
1834	!
1835	!-- If 10-m level is below the first grid level, MOST is
1836	!-- used for calculation of 10-m temperature.
1837	IF ( surf%z_mo(m) > 10.0_wp ) THEN
1838	uv_10m(j,i) = surf%us(m) / kappa &
1839	* ( LOG( 10.0_wp / surf%z0(m) ) &
1840	- psi_m( 10.0_wp / surf%ol(m) ) &
1841	+ psi_m( surf%z0(m) / surf%ol(m) ) )
1842	!
1843	!-- If 10-m level is above the first grid level, 10-m wind speed
1844	!-- is linearly interpolated between the two nearest vertical grid
1845	!-- levels. Note, since 10-m temperature is only computed for
1846	!-- horizontal upward-facing surfaces, only a vertical
1847	!-- interpolation is necessary.
1848	ELSE
1849	!
1850	!-- zw(k-1) defines the height of the surface.
1851	kk = k
1852	DO WHILE ( zu(kk) - zw(k-1) < 10.0_wp .AND. kk <= nzt )
1853	kk = kk + 1
1854	ENDDO
1855	!
1856	!-- kk defines the index of the first grid level >= 10m.
1857	uv_l = SQRT( ( 0.5_wp * ( u(kk-1,j,i) + u(kk-1,j,i+1) ) )**2 &
1858	+ ( 0.5_wp * ( v(kk-1,j,i) + v(kk-1,j+1,i) ) )**2 )
1859
1860	uv_u = SQRT( ( 0.5_wp * ( u(kk,j,i) + u(kk,j,i+1) ) )**2 &
1861	+ ( 0.5_wp * ( v(kk,j,i) + v(kk,j+1,i) ) )**2 )
1862
1863	uv_10m(j,i) = uv_l + ( zw(k-1) + 10.0_wp - zu(kk-1) ) * &
1864	( uv_u - uv_l ) / &
1865	( zu(kk) - zu(kk-1) )
1866
1867	ENDIF
1868
1869	ENDDO
1870
1871	END SUBROUTINE calc_wind_10m
1872
1873	END SUBROUTINE doq_calculate
1874
1875
1876	!------------------------------------------------------------------------------!
1877	! Description:
1878	! ------------
1879	!> Preparation of the diagnostic output, counting of the module-specific
1880	!> output quantities and gathering of the output names.
1881	!------------------------------------------------------------------------------!
1882	SUBROUTINE doq_prepare
1883
1884	USE control_parameters, &
1885	ONLY: do2d, do3d, domask, masks
1886
1887	IMPLICIT NONE
1888
1889	CHARACTER (LEN=varnamelength), DIMENSION(0:1,500) :: do2d_var = ' ' !<
1890	!< label array for 2d output quantities
1891
1892	INTEGER(iwp) :: av !< index defining type of output, av=0 instantaneous, av=1 averaged
1893	INTEGER(iwp) :: ivar !< loop index
1894	INTEGER(iwp) :: ivar_all !< loop index
1895	INTEGER(iwp) :: l !< index for cutting string
1896	INTEGER(iwp) :: mid !< masked output running index
1897
1898	prepared_diagnostic_output_quantities = .FALSE.
1899
1900	ivar = 1
1901	ivar_all = 1
1902
1903	DO av = 0, 1
1904	!
1905	!-- Remove _xy, _xz, or _yz from string
1906	l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
1907	do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
1908	!
1909	!-- Gather 2d output quantity names.
1910	!-- Check for double occurrence of output quantity, e.g. by _xy,
1911	!-- _yz, _xz.
1912	DO WHILE ( do2d_var(av,ivar)(1:1) /= ' ' )
1913	IF ( .NOT. ANY( do_all == do2d_var(av,ivar) ) ) THEN
1914	do_all(ivar_all) = do2d_var(av,ivar)
1915	ENDIF
1916	ivar = ivar + 1
1917	ivar_all = ivar_all + 1
1918	l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
1919	do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
1920	ENDDO
1921
1922	ivar = 1
1923	!
1924	!-- Gather 3d output quantity names
1925	DO WHILE ( do3d(av,ivar)(1:1) /= ' ' )
1926	do_all(ivar_all) = do3d(av,ivar)
1927	ivar = ivar + 1
1928	ivar_all = ivar_all + 1
1929	ENDDO
1930
1931	ivar = 1
1932	!
1933	!-- Gather masked output quantity names. Also check for double output
1934	!-- e.g. by different masks.
1935	DO mid = 1, masks
1936	DO WHILE ( domask(mid,av,ivar)(1:1) /= ' ' )
1937	IF ( .NOT. ANY( do_all == domask(mid,av,ivar) ) ) THEN
1938	do_all(ivar_all) = domask(mid,av,ivar)
1939	ENDIF
1940
1941	ivar = ivar + 1
1942	ivar_all = ivar_all + 1
1943	ENDDO
1944	ivar = 1
1945	ENDDO
1946
1947	ENDDO
1948
1949	prepared_diagnostic_output_quantities = .TRUE.
1950
1951	END SUBROUTINE doq_prepare
1952
1953	!------------------------------------------------------------------------------!
1954	! Description:
1955	! ------------
1956	!> Subroutine reads local (subdomain) restart data
1957	!------------------------------------------------------------------------------!
1958	SUBROUTINE doq_rrd_local_ftn( k, nxlf, nxlc, nxl_on_file, nxrf, nxrc, &
1959	nxr_on_file, nynf, nync, nyn_on_file, nysf, &
1960	nysc, nys_on_file, tmp_2d, tmp_3d, found )
1961
1962	USE control_parameters
1963
1964	USE indices
1965
1966	USE kinds
1967
1968	USE pegrid
1969
1970
1971	IMPLICIT NONE
1972
1973	INTEGER(iwp) :: k !<
1974	INTEGER(iwp) :: nxlc !<
1975	INTEGER(iwp) :: nxlf !<
1976	INTEGER(iwp) :: nxl_on_file !<
1977	INTEGER(iwp) :: nxrc !<
1978	INTEGER(iwp) :: nxrf !<
1979	INTEGER(iwp) :: nxr_on_file !<
1980	INTEGER(iwp) :: nync !<
1981	INTEGER(iwp) :: nynf !<
1982	INTEGER(iwp) :: nyn_on_file !<
1983	INTEGER(iwp) :: nysc !<
1984	INTEGER(iwp) :: nysf !<
1985	INTEGER(iwp) :: nys_on_file !<
1986
1987	LOGICAL, INTENT(OUT) :: found
1988
1989	REAL(wp), DIMENSION(nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_2d !<
1990
1991	REAL(wp), DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d !<
1992
1993
1994	found = .TRUE.
1995
1996	SELECT CASE ( restart_string(1:length) )
1997
1998	CASE ( 'pt_2m_av' )
1999	IF ( .NOT. ALLOCATED( pt_2m_av ) ) THEN
2000	ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
2001	ENDIF
2002	IF ( k == 1 ) READ ( 13 ) tmp_2d
2003	pt_2m_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2004	tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2005
2006	CASE ( 'ti_av' )
2007	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
2008	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2009	ENDIF
2010	IF ( k == 1 ) READ ( 13 ) tmp_3d
2011	ti_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2012	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2013
2014	CASE ( 'u_center_av' )
2015	IF ( .NOT. ALLOCATED( u_center_av ) ) THEN
2016	ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2017	ENDIF
2018	IF ( k == 1 ) READ ( 13 ) tmp_3d
2019	u_center_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2020	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2021
2022	CASE ( 'uu_av' )
2023	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
2024	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2025	ENDIF
2026	IF ( k == 1 ) READ ( 13 ) tmp_3d
2027	uu_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2028	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2029
2030	CASE ( 'uv_10m_av' )
2031	IF ( .NOT. ALLOCATED( uv_10m_av ) ) THEN
2032	ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
2033	ENDIF
2034	IF ( k == 1 ) READ ( 13 ) tmp_2d
2035	uv_10m_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2036	tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2037
2038	CASE ( 'v_center_av' )
2039	IF ( .NOT. ALLOCATED( v_center_av ) ) THEN
2040	ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2041	ENDIF
2042	IF ( k == 1 ) READ ( 13 ) tmp_3d
2043	v_center_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2044	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2045
2046	CASE ( 'vv_av' )
2047	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
2048	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2049	ENDIF
2050	IF ( k == 1 ) READ ( 13 ) tmp_3d
2051	vv_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2052	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2053
2054	CASE ( 'wtheta_av' )
2055	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
2056	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2057	ENDIF
2058	IF ( k == 1 ) READ ( 13 ) tmp_3d
2059	wtheta_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2060	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2061
2062	CASE ( 'wq_av' )
2063	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
2064	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2065	ENDIF
2066	IF ( k == 1 ) READ ( 13 ) tmp_3d
2067	wq_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2068	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2069
2070	CASE ( 'wspeed_av' )
2071	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
2072	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2073	ENDIF
2074	IF ( k == 1 ) READ ( 13 ) tmp_3d
2075	wspeed_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2076	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2077
2078	CASE ( 'wu_av' )
2079	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
2080	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2081	ENDIF
2082	IF ( k == 1 ) READ ( 13 ) tmp_3d
2083	wu_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2084	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2085
2086	CASE ( 'wv_av' )
2087	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
2088	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2089	ENDIF
2090	IF ( k == 1 ) READ ( 13 ) tmp_3d
2091	wv_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2092	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2093
2094	CASE ( 'ww_av' )
2095	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
2096	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2097	ENDIF
2098	IF ( k == 1 ) READ ( 13 ) tmp_3d
2099	ww_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
2100	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2101
2102
2103	CASE DEFAULT
2104
2105	found = .FALSE.
2106
2107	END SELECT
2108
2109	END SUBROUTINE doq_rrd_local_ftn
2110
2111	!------------------------------------------------------------------------------!
2112	! Description:
2113	! ------------
2114	!> Read module-specific local restart data arrays (MPI-IO).
2115	!------------------------------------------------------------------------------!
2116	SUBROUTINE doq_rrd_local_mpi
2117
2118	LOGICAL :: array_found !< control flad indicating whether the array is found in the file or not
2119
2120
2121	CALL rd_mpi_io_check_array( 'pt_2m_av' , found = array_found )
2122	IF ( array_found ) THEN
2123	IF ( .NOT. ALLOCATED( pt_2m_av ) ) ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
2124	CALL rrd_mpi_io( 'pt_2m_av', pt_2m_av )
2125	ENDIF
2126
2127	CALL rd_mpi_io_check_array( 'ti_av' , found = array_found )
2128	IF ( array_found ) THEN
2129	IF ( .NOT. ALLOCATED( ti_av ) ) ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2130	CALL rrd_mpi_io( 'ti_av', ti_av )
2131	ENDIF
2132
2133	CALL rd_mpi_io_check_array( 'u_center_av' , found = array_found )
2134	IF ( array_found ) THEN
2135	IF ( .NOT. ALLOCATED( u_center_av ) ) ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2136	CALL rrd_mpi_io( 'u_center_av', u_center_av )
2137	ENDIF
2138
2139	CALL rd_mpi_io_check_array( 'uu_av' , found = array_found )
2140	IF ( array_found ) THEN
2141	IF ( .NOT. ALLOCATED( uu_av ) ) ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2142	CALL rrd_mpi_io( 'uu_av', uu_av )
2143	ENDIF
2144
2145	CALL rd_mpi_io_check_array( 'uv_10m_av' , found = array_found )
2146	IF ( array_found ) THEN
2147	IF ( .NOT. ALLOCATED( uv_10m_av ) ) ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
2148	CALL rrd_mpi_io( 'uv_10m_av', uv_10m_av )
2149	ENDIF
2150
2151	CALL rd_mpi_io_check_array( 'v_center_av' , found = array_found )
2152	IF ( array_found ) THEN
2153	IF ( .NOT. ALLOCATED( v_center_av ) ) ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2154	CALL rrd_mpi_io( 'v_center_av', v_center_av )
2155	ENDIF
2156
2157	CALL rd_mpi_io_check_array( 'vv_av' , found = array_found )
2158	IF ( array_found ) THEN
2159	IF ( .NOT. ALLOCATED( vv_av ) ) ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2160	CALL rrd_mpi_io( 'vv_av', vv_av )
2161	ENDIF
2162
2163	CALL rd_mpi_io_check_array( 'ww_av' , found = array_found )
2164	IF ( array_found ) THEN
2165	IF ( .NOT. ALLOCATED( ww_av ) ) ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2166	CALL rrd_mpi_io( 'ww_av', ww_av )
2167	ENDIF
2168
2169	CALL rd_mpi_io_check_array( 'wu_av' , found = array_found )
2170	IF ( array_found ) THEN
2171	IF ( .NOT. ALLOCATED( wu_av ) ) ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2172	CALL rrd_mpi_io( 'wu_av', wu_av )
2173	ENDIF
2174
2175	CALL rd_mpi_io_check_array( 'wv_av' , found = array_found )
2176	IF ( array_found ) THEN
2177	IF ( .NOT. ALLOCATED( wv_av ) ) ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2178	CALL rrd_mpi_io( 'wv_av', wv_av )
2179	ENDIF
2180
2181	CALL rd_mpi_io_check_array( 'wtheta_av' , found = array_found )
2182	IF ( array_found ) THEN
2183	IF ( .NOT. ALLOCATED( wtheta_av ) ) ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2184	CALL rrd_mpi_io( 'wtheta_av', wtheta_av )
2185	ENDIF
2186
2187	CALL rd_mpi_io_check_array( 'wq_av' , found = array_found )
2188	IF ( array_found ) THEN
2189	IF ( .NOT. ALLOCATED( wq_av ) ) ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2190	CALL rrd_mpi_io( 'wq_av', wq_av )
2191	ENDIF
2192
2193	CALL rd_mpi_io_check_array( 'wspeed_av' , found = array_found )
2194	IF ( array_found ) THEN
2195	IF ( .NOT. ALLOCATED( wspeed_av ) ) ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2196	CALL rrd_mpi_io( 'wspeed_av', wspeed_av )
2197	ENDIF
2198
2199	END SUBROUTINE doq_rrd_local_mpi
2200
2201	!------------------------------------------------------------------------------!
2202	! Description:
2203	! ------------
2204	!> Subroutine writes local (subdomain) restart data
2205	!------------------------------------------------------------------------------!
2206	SUBROUTINE doq_wrd_local
2207
2208
2209	IMPLICIT NONE
2210
2211	IF ( TRIM( restart_data_format_output ) == 'fortran_binary' ) THEN
2212
2213	IF ( ALLOCATED( pt_2m_av ) ) THEN
2214	CALL wrd_write_string( 'pt_2m_av' )
2215	WRITE ( 14 ) pt_2m_av
2216	ENDIF
2217
2218	IF ( ALLOCATED( ti_av ) ) THEN
2219	CALL wrd_write_string( 'ti_av' )
2220	WRITE ( 14 ) ti_av
2221	ENDIF
2222
2223	IF ( ALLOCATED( u_center_av ) ) THEN
2224	CALL wrd_write_string( 'u_center_av' )
2225	WRITE ( 14 ) u_center_av
2226	ENDIF
2227
2228	IF ( ALLOCATED( uu_av ) ) THEN
2229	CALL wrd_write_string( 'uu_av' )
2230	WRITE ( 14 ) uu_av
2231	ENDIF
2232
2233	IF ( ALLOCATED( uv_10m_av ) ) THEN
2234	CALL wrd_write_string( 'uv_10m_av' )
2235	WRITE ( 14 ) uv_10m_av
2236	ENDIF
2237
2238	IF ( ALLOCATED( v_center_av ) ) THEN
2239	CALL wrd_write_string( 'v_center_av' )
2240	WRITE ( 14 ) v_center_av
2241	ENDIF
2242
2243	IF ( ALLOCATED( vv_av ) ) THEN
2244	CALL wrd_write_string( 'vv_av' )
2245	WRITE ( 14 ) vv_av
2246	ENDIF
2247
2248	IF ( ALLOCATED( ww_av ) ) THEN
2249	CALL wrd_write_string( 'ww_av' )
2250	WRITE ( 14 ) ww_av
2251	ENDIF
2252
2253	IF ( ALLOCATED( wu_av ) ) THEN
2254	CALL wrd_write_string( 'wu_av' )
2255	WRITE ( 14 ) wu_av
2256	ENDIF
2257
2258	IF ( ALLOCATED( wv_av ) ) THEN
2259	CALL wrd_write_string( 'wv_av' )
2260	WRITE ( 14 ) wv_av
2261	ENDIF
2262
2263	IF ( ALLOCATED( wtheta_av ) ) THEN
2264	CALL wrd_write_string( 'wtheta_av' )
2265	WRITE ( 14 ) wtheta_av
2266	ENDIF
2267
2268	IF ( ALLOCATED( wq_av ) ) THEN
2269	CALL wrd_write_string( 'wq_av' )
2270	WRITE ( 14 ) wq_av
2271	ENDIF
2272
2273	IF ( ALLOCATED( wspeed_av ) ) THEN
2274	CALL wrd_write_string( 'wspeed_av' )
2275	WRITE ( 14 ) wspeed_av
2276	ENDIF
2277
2278	ELSEIF ( restart_data_format_output(1:3) == 'mpi' ) THEN
2279
2280	IF ( ALLOCATED( pt_2m_av ) ) CALL wrd_mpi_io( 'pt_2m_av', pt_2m_av )
2281	IF ( ALLOCATED( ti_av ) ) CALL wrd_mpi_io( 'ti_av', ti_av )
2282	IF ( ALLOCATED( u_center_av ) ) CALL wrd_mpi_io( 'u_center_av', u_center_av )
2283	IF ( ALLOCATED( uu_av ) ) CALL wrd_mpi_io( 'uu_av', uu_av )
2284	IF ( ALLOCATED( uv_10m_av ) ) CALL wrd_mpi_io( 'uv_10m_av', uv_10m_av )
2285	IF ( ALLOCATED( vv_av ) ) CALL wrd_mpi_io( 'vv_av', vv_av )
2286	IF ( ALLOCATED( v_center_av ) ) CALL wrd_mpi_io( 'v_center_av', v_center_av )
2287	IF ( ALLOCATED( wtheta_av ) ) CALL wrd_mpi_io( 'wtheta_av', wtheta_av )
2288	IF ( ALLOCATED( wq_av ) ) CALL wrd_mpi_io( 'wq_av', wq_av )
2289	IF ( ALLOCATED( wspeed_av ) ) CALL wrd_mpi_io( 'wspeed_av', wspeed_av )
2290	IF ( ALLOCATED( wu_av ) ) CALL wrd_mpi_io( 'wu_av', wu_av )
2291	IF ( ALLOCATED( wv_av ) ) CALL wrd_mpi_io( 'wv_av', wv_av )
2292	IF ( ALLOCATED( ww_av ) ) CALL wrd_mpi_io( 'ww_av', ww_av )
2293
2294	ENDIF
2295
2296	END SUBROUTINE doq_wrd_local
2297
2298
2299	END MODULE diagnostic_output_quantities_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |