Home

Context Navigation

source: palm/trunk/SOURCE/data_output_dvrp.f90 @ 2309

Last change on this file since 2309 was 2274, checked in by Giersch, 7 years ago
Changed error messages
Property svn:keywords set to `Id`
File size: 20.7 KB

Line
1	!> @file data_output_dvrp.f90
2	!------------------------------------------------------------------------------!
3	! This file is part of PALM.
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-2017 Leibniz Universitaet Hannover
18	!------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! -----------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: data_output_dvrp.f90 2274 2017-06-09 13:27:48Z gronemeier $
27	! Changed error messages
28	!
29	! 2101 2017-01-05 16:42:31Z suehring
30	!
31	! 2000 2016-08-20 18:09:15Z knoop
32	! Forced header and separation lines into 80 columns
33	!
34	! 1960 2016-07-12 16:34:24Z suehring
35	! Separate humidity and passive scalar
36	!
37	! 1873 2016-04-18 14:50:06Z maronga
38	! Module renamed (removed _mod)
39	!
40	!
41	! 1850 2016-04-08 13:29:27Z maronga
42	! Module renamed
43	!
44	!
45	! 1822 2016-04-07 07:49:42Z hoffmann
46	! Particles and tails removed.
47	!
48	! 1682 2015-10-07 23:56:08Z knoop
49	! Code annotations made doxygen readable
50	!
51	! 1353 2014-04-08 15:21:23Z heinze
52	! REAL constants provided with KIND-attribute
53	!
54	! 1346 2014-03-27 13:18:20Z heinze
55	! Bugfix: REAL constants provided with KIND-attribute especially in call of
56	! intrinsic function like MAX, MIN, SIGN
57	!
58	! 1320 2014-03-20 08:40:49Z raasch
59	! ONLY-attribute added to USE-statements,
60	! kind-parameters added to all INTEGER and REAL declaration statements,
61	! kinds are defined in new module kinds,
62	! revision history before 2012 removed,
63	! comment fields (!:) to be used for variable explanations added to
64	! all variable declaration statements
65	!
66	! 1318 2014-03-17 13:35:16Z raasch
67	! module interfaces removed
68	!
69	! 1036 2012-10-22 13:43:42Z raasch
70	! code put under GPL (PALM 3.9)
71	!
72	! 828 2012-02-21 12:00:36Z raasch
73	! particle feature color renamed class
74	!
75	! Revision 1.1 2000/04/27 06:27:17 raasch
76	! Initial revision
77	!
78	!
79	! Description:
80	! ------------
81	!> Plot of isosurface and slicers with dvrp-software
82	!------------------------------------------------------------------------------!
83	MODULE dvrp_color
84
85
86	USE dvrp_variables
87
88	USE kinds
89
90	IMPLICIT NONE
91
92	CONTAINS
93
94	!------------------------------------------------------------------------------!
95	! Description:
96	! ------------
97	!> @todo Missing subroutine description.
98	!------------------------------------------------------------------------------!
99	SUBROUTINE color_dvrp( value, color )
100
101	REAL(wp), INTENT(IN) :: value !<
102	REAL(wp), INTENT(OUT) :: color(4) !<
103
104	REAL(wp) :: scale !<
105
106	scale = ( value - slicer_range_limits_dvrp(1,islice_dvrp) ) / &
107	( slicer_range_limits_dvrp(2,islice_dvrp) - &
108	slicer_range_limits_dvrp(1,islice_dvrp) )
109
110	scale = MODULO( 180.0_wp + 180.0_wp * scale, 360.0_wp )
111
112	color = (/ scale, 0.5_wp, 1.0_wp, 0.0_wp /)
113
114	END SUBROUTINE color_dvrp
115
116	END MODULE dvrp_color
117
118
119	RECURSIVE SUBROUTINE data_output_dvrp
120
121	#if defined( __dvrp_graphics )
122
123	USE arrays_3d, &
124	ONLY: p, pt, q, ql, s, ts, u, us, v, w, zu
125
126	USE cloud_parameters, &
127	ONLY: l_d_cp, pt_d_t
128
129	USE constants, &
130	ONLY: pi
131
132	USE control_parameters, &
133	ONLY: cloud_droplets, cloud_physics, do2d, do3d, humidity, ibc_uv_b, &
134	message_string, nz_do3d, passive_scalar, simulated_time, &
135	threshold
136
137	USE cpulog, &
138	ONLY: log_point, log_point_s, cpu_log
139
140	USE DVRP
141
142	USE dvrp_color
143
144	USE dvrp_variables
145
146	USE grid_variables, &
147	ONLY: dx, dy
148
149	USE indices, &
150	ONLY: nxl, nxr, nyn, nys, nzb
151
152	USE kinds
153
154	USE pegrid
155
156	IMPLICIT NONE
157
158	CHARACTER (LEN=2) :: section_chr !<
159	CHARACTER (LEN=6) :: output_variable !<
160
161	INTEGER(iwp) :: c_mode !<
162	INTEGER(iwp) :: c_size_x !<
163	INTEGER(iwp) :: c_size_y !<
164	INTEGER(iwp) :: c_size_z !<
165	INTEGER(iwp) :: dvrp_nop !<
166	INTEGER(iwp) :: dvrp_not !<
167	INTEGER(iwp) :: gradient_normals !<
168	INTEGER(iwp) :: i !<
169	INTEGER(iwp) :: ip !<
170	INTEGER(iwp) :: j !<
171	INTEGER(iwp) :: jp !<
172	INTEGER(iwp) :: k !<
173	INTEGER(iwp) :: l !<
174	INTEGER(iwp) :: m !<
175	INTEGER(iwp) :: n !<
176	INTEGER(iwp) :: n_isosurface !<
177	INTEGER(iwp) :: n_slicer !<
178	INTEGER(iwp) :: nn !<
179	INTEGER(iwp) :: section_mode !<
180	INTEGER(iwp) :: vn !<
181	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: p_c !<
182	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: p_t !<
183
184	LOGICAL, DIMENSION(:), ALLOCATABLE :: dvrp_mask !<
185
186	REAL(sp) :: slicer_position !<
187	REAL(sp) :: tmp_alpha !<
188	REAL(sp) :: tmp_alpha_w !<
189	REAL(sp) :: tmp_b !<
190	REAL(sp) :: tmp_c_alpha !<
191	REAL(sp) :: tmp_g !<
192	REAL(sp) :: tmp_norm !<
193	REAL(sp) :: tmp_pos !<
194	REAL(sp) :: tmp_r !<
195	REAL(sp) :: tmp_t !<
196	REAL(sp) :: tmp_th !<
197	REAL(sp), DIMENSION(:), ALLOCATABLE :: psize !<
198	REAL(sp), DIMENSION(:), ALLOCATABLE :: p_x !<
199	REAL(sp), DIMENSION(:), ALLOCATABLE :: p_y !<
200	REAL(sp), DIMENSION(:), ALLOCATABLE :: p_z !<
201	REAL(sp), DIMENSION(:,:,:), ALLOCATABLE :: local_pf !<
202	REAL(sp), DIMENSION(:,:,:,:), ALLOCATABLE :: local_pfi !<
203
204
205	CALL cpu_log( log_point(27), 'data_output_dvrp', 'start' )
206
207	!
208	!-- Loop over all output modes choosed
209	m = 1
210	n_isosurface = 0 ! isosurface counter (for threshold values and color)
211	n_slicer = 0 ! slice plane counter (for range of values)
212	DO WHILE ( mode_dvrp(m) /= ' ' )
213	!
214	!-- Update of the steering variables
215	IF ( .NOT. lock_steering_update ) THEN
216	!
217	!-- Set lock to avoid recursive calls of DVRP_STEERING_UPDATE
218	lock_steering_update = .TRUE.
219	! CALL DVRP_STEERING_UPDATE( m-1, data_output_dvrp )
220	lock_steering_update = .FALSE.
221	ENDIF
222
223	!
224	!-- Determine the variable which shall be plotted (in case of slicers or
225	!-- isosurfaces)
226	IF ( mode_dvrp(m)(1:10) == 'isosurface' ) THEN
227	READ ( mode_dvrp(m), '(10X,I2)' ) vn
228	output_variable = do3d(0,vn)
229	n_isosurface = n_isosurface + 1
230	ELSEIF ( mode_dvrp(m)(1:6) == 'slicer' ) THEN
231	READ ( mode_dvrp(m), '(6X,I2)' ) vn
232	output_variable = do2d(0,vn)
233	l = MAX( 2, LEN_TRIM( do2d(0,vn) ) )
234	section_chr = do2d(0,vn)(l-1:l)
235	SELECT CASE ( section_chr )
236	CASE ( 'xy' )
237	section_mode = 2
238	slicer_position = zu(MIN( slicer_position_dvrp(m), nz_do3d ))
239	CASE ( 'xz' )
240	section_mode = 1
241	slicer_position = slicer_position_dvrp(m) * dy
242	CASE ( 'yz' )
243	section_mode = 0
244	slicer_position = slicer_position_dvrp(m) * dx
245	END SELECT
246	ENDIF
247
248	!
249	!-- Select the plot mode (in case of isosurface or slicer only if user has
250	!-- defined a variable which shall be plotted; otherwise do nothing)
251	IF ( ( mode_dvrp(m)(1:10) == 'isosurface' .OR. &
252	mode_dvrp(m)(1:6) == 'slicer' ) &
253	.AND. output_variable /= ' ' ) THEN
254
255	!
256	!-- Create an intermediate array, properly dimensioned for plot-output
257	ALLOCATE( local_pf(nxl_dvrp:nxr_dvrp+1,nys_dvrp:nyn_dvrp+1, &
258	nzb:nz_do3d) )
259
260	!
261	!-- Move original array to intermediate array
262	IF ( dvrp_overlap ) THEN
263
264	SELECT CASE ( output_variable )
265
266	CASE ( 'u', 'u_xy', 'u_xz', 'u_yz' )
267	DO i = nxl_dvrp, nxr_dvrp+1
268	DO j = nys_dvrp, nyn_dvrp+1
269	DO k = nzb, nz_do3d
270	local_pf(i,j,k) = u(k,j,i)
271	ENDDO
272	ENDDO
273	ENDDO
274	!
275	!-- Replace mirrored values at lower surface by real surface
276	!-- values
277	IF ( output_variable == 'u_xz' .OR. &
278	output_variable == 'u_yz' ) THEN
279	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0_wp
280	ENDIF
281
282
283	CASE ( 'v', 'v_xy', 'v_xz', 'v_yz' )
284	DO i = nxl_dvrp, nxr_dvrp+1
285	DO j = nys_dvrp, nyn_dvrp+1
286	DO k = nzb, nz_do3d
287	local_pf(i,j,k) = v(k,j,i)
288	ENDDO
289	ENDDO
290	ENDDO
291	!
292	!-- Replace mirrored values at lower surface by real surface
293	!-- values
294	IF ( output_variable == 'v_xz' .OR. &
295	output_variable == 'v_yz' ) THEN
296	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0_wp
297	ENDIF
298
299	CASE ( 'w', 'w_xy', 'w_xz', 'w_yz' )
300	DO i = nxl_dvrp, nxr_dvrp+1
301	DO j = nys_dvrp, nyn_dvrp+1
302	DO k = nzb, nz_do3d
303	local_pf(i,j,k) = w(k,j,i)
304	ENDDO
305	ENDDO
306	ENDDO
307	! Averaging for Langmuir circulation
308	! DO k = nzb, nz_do3d
309	! DO j = nys_dvrp+1, nyn_dvrp
310	! DO i = nxl_dvrp, nxr_dvrp+1
311	! local_pf(i,j,k) = 0.25 * local_pf(i,j-1,k) + &
312	! 0.50 * local_pf(i,j,k) + &
313	! 0.25 * local_pf(i,j+1,k)
314	! ENDDO
315	! ENDDO
316	! ENDDO
317
318	CASE ( 'p', 'p_xy', 'p_xz', 'p_yz' )
319	DO i = nxl_dvrp, nxr_dvrp+1
320	DO j = nys_dvrp, nyn_dvrp+1
321	DO k = nzb, nz_do3d
322	local_pf(i,j,k) = p(k,j,i)
323	ENDDO
324	ENDDO
325	ENDDO
326
327	CASE ( 'pt', 'pt_xy', 'pt_xz', 'pt_yz' )
328	IF ( .NOT. cloud_physics ) THEN
329	DO i = nxl_dvrp, nxr_dvrp+1
330	DO j = nys_dvrp, nyn_dvrp+1
331	DO k = nzb, nz_do3d
332	local_pf(i,j,k) = pt(k,j,i)
333	ENDDO
334	ENDDO
335	ENDDO
336	ELSE
337	DO i = nxl_dvrp, nxr_dvrp+1
338	DO j = nys_dvrp, nyn_dvrp+1
339	DO k = nzb, nz_do3d
340	local_pf(i,j,k) = pt(k,j,i) + l_d_cp * &
341	pt_d_t(k) * ql(k,j,i)
342	ENDDO
343	ENDDO
344	ENDDO
345	ENDIF
346
347	CASE ( 'q', 'q_xy', 'q_xz', 'q_yz' )
348	IF ( humidity ) THEN
349	DO i = nxl_dvrp, nxr_dvrp+1
350	DO j = nys_dvrp, nyn_dvrp+1
351	DO k = nzb, nz_do3d
352	local_pf(i,j,k) = q(k,j,i)
353	ENDDO
354	ENDDO
355	ENDDO
356	ELSE
357	message_string = 'if humidity = ' // &
358	'.FALSE. output of ' // TRIM( output_variable ) // &
359	'is not provided'
360	CALL message( 'data_output_dvrp', 'PA0183',&
361	0, 0, 0, 6, 0 )
362	ENDIF
363
364	CASE ( 'ql', 'ql_xy', 'ql_xz', 'ql_yz' )
365	IF ( cloud_physics .OR. cloud_droplets ) THEN
366	DO i = nxl_dvrp, nxr_dvrp+1
367	DO j = nys_dvrp, nyn_dvrp+1
368	DO k = nzb, nz_do3d
369	local_pf(i,j,k) = ql(k,j,i)
370	ENDDO
371	ENDDO
372	ENDDO
373	ELSE
374	message_string = 'if cloud_physics = .FALSE. and ' // &
375	'cloud_droplets = .FALSE. '
376	'output of ' // TRIM( output_variable) // &
377	'is not provided'
378	CALL message( 'data_output_dvrp', 'PA0184',&
379	0, 0, 0, 6, 0 )
380	ENDIF
381
382	CASE ( 's', 's_xy', 's_xz', 's_yz' )
383	IF ( passive_scalar ) THEN
384	DO i = nxl_dvrp, nxr_dvrp+1
385	DO j = nys_dvrp, nyn_dvrp+1
386	DO k = nzb, nz_do3d
387	local_pf(i,j,k) = s(k,j,i)
388	ENDDO
389	ENDDO
390	ENDDO
391	ELSE
392	message_string = 'if passive_scalar = ' // &
393	'.FALSE. output of ' // TRIM( output_variable ) // &
394	'is not provided'
395	CALL message( 'data_output_dvrp', 'PA0183',&
396	0, 0, 0, 6, 0 )
397	ENDIF
398
399	CASE ( 'u*_xy' )
400	DO i = nxl_dvrp, nxr_dvrp+1
401	DO j = nys_dvrp, nyn_dvrp+1
402	local_pf(i,j,nzb+1) = us(j,i)
403	ENDDO
404	ENDDO
405	slicer_position = zu(nzb+1)
406
407	CASE ( 't*_xy' )
408	DO i = nxl_dvrp, nxr_dvrp+1
409	DO j = nys_dvrp, nyn_dvrp+1
410	local_pf(i,j,nzb+1) = ts(j,i)
411	ENDDO
412	ENDDO
413	slicer_position = zu(nzb+1)
414
415
416	CASE DEFAULT
417	!
418	!-- The DEFAULT case is reached either if output_variable
419	!-- contains unsupported variable or if the user has coded a
420	!-- special case in the user interface. There, the subroutine
421	!-- user_data_output_dvrp checks which of these two conditions
422	!-- applies.
423	CALL user_data_output_dvrp( output_variable, local_pf )
424
425
426	END SELECT
427
428	ELSE
429	!
430	!-- No overlap of clipping domain with the current subdomain
431	DO i = nxl_dvrp, nxr_dvrp+1
432	DO j = nys_dvrp, nyn_dvrp+1
433	DO k = nzb, nz_do3d
434	local_pf(i,j,k) = 0.0_wp
435	ENDDO
436	ENDDO
437	ENDDO
438
439	ENDIF
440
441	IF ( mode_dvrp(m)(1:10) == 'isosurface' ) THEN
442
443	!
444	!-- DVRP-Calls for plotting isosurfaces:
445	CALL cpu_log( log_point_s(26), 'dvrp_isosurface', 'start' )
446
447	!
448	!-- Definition of isosurface color
449	tmp_r = isosurface_color(1,n_isosurface)
450	tmp_g = isosurface_color(2,n_isosurface)
451	tmp_b = isosurface_color(3,n_isosurface)
452	tmp_t = 0.0_wp
453	CALL DVRP_MATERIAL_RGB( m-1, 1, tmp_r, tmp_g, tmp_b, tmp_t )
454
455	!
456	!-- Compute and plot isosurface in dvr-format
457	CALL DVRP_DATA( m-1, local_pf, 1, nx_dvrp, ny_dvrp, nz_dvrp, &
458	cyclic_dvrp, cyclic_dvrp, cyclic_dvrp )
459
460	c_size_x = vc_size_x; c_size_y = vc_size_y; c_size_z = vc_size_z
461	CALL DVRP_CLUSTER_SIZE( m-1, c_size_x, c_size_y, c_size_z )
462
463	c_mode = vc_mode
464	CALL DVRP_CLUSTERING_MODE( m-1, c_mode )
465
466	gradient_normals = vc_gradient_normals
467	CALL DVRP_GRADIENTNORMALS( m-1, gradient_normals )
468
469	!
470	!-- A seperate procedure for setting vc_alpha will be in the next
471	!-- version of libDVRP
472	tmp_c_alpha = vc_alpha
473	CALL DVRP_THRESHOLD( -(m-1)-1, tmp_c_alpha )
474
475	IF ( dvrp_overlap ) THEN
476	tmp_th = threshold(n_isosurface)
477	ELSE
478	tmp_th = 1.0_wp ! nothing is plotted because array values are 0
479	ENDIF
480
481	CALL DVRP_THRESHOLD( m-1, tmp_th )
482
483	CALL DVRP_VISUALIZE( m-1, 21, dvrp_filecount )
484
485	CALL cpu_log( log_point_s(26), 'dvrp_isosurface', 'stop' )
486
487	ELSEIF ( mode_dvrp(m)(1:6) == 'slicer' ) THEN
488
489	!
490	!-- DVRP-Calls for plotting slicers:
491	CALL cpu_log( log_point_s(27), 'dvrp_slicer', 'start' )
492
493	!
494	!-- Material and color definitions
495	tmp_r = 0.0_wp; tmp_g = 0.0_wp; tmp_b = 0.0_wp; tmp_t = 0.0_wp
496	CALL DVRP_MATERIAL_RGB( m-1, 1, tmp_r, tmp_g, tmp_b, tmp_t )
497
498	n_slicer = n_slicer + 1
499
500	!
501	!-- Using dolorfunction has not been properly tested
502	! islice_dvrp = n_slicer
503	! CALL DVRP_COLORFUNCTION( m-1, DVRP_CM_HLS, 25, &
504	! slicer_range_limits_dvrp(:,n_slicer), &
505	! color_dvrp )
506
507	!
508	!-- Set interval of values defining the colortable
509	CALL set_slicer_attributes_dvrp( n_slicer )
510
511	!
512	!-- Create user-defined colortable
513	CALL user_dvrp_coltab( 'slicer', output_variable )
514
515	CALL DVRP_COLORTABLE_HLS( m-1, 1, interval_values_dvrp, &
516	interval_h_dvrp, interval_l_dvrp, &
517	interval_s_dvrp, interval_a_dvrp )
518
519	!
520	!-- Compute and plot slicer in dvr-format
521	CALL DVRP_DATA( m-1, local_pf, 1, nx_dvrp, ny_dvrp, nz_dvrp, &
522	cyclic_dvrp, cyclic_dvrp, cyclic_dvrp )
523	tmp_pos = slicer_position
524	CALL DVRP_SLICER( m-1, section_mode, tmp_pos )
525
526	CALL DVRP_VISUALIZE( m-1, 2, dvrp_filecount )
527
528	CALL cpu_log( log_point_s(27), 'dvrp_slicer', 'stop' )
529
530	ENDIF
531
532	DEALLOCATE( local_pf )
533
534	ELSEIF ( mode_dvrp(m)(1:9) == 'pathlines' ) THEN
535
536	ALLOCATE( local_pfi(4,nxl:nxr+1,nys:nyn+1,nzb:nz_do3d) )
537	DO i = nxl, nxr+1
538	DO j = nys, nyn+1
539	DO k = nzb, nz_do3d
540	local_pfi(1,i,j,k) = u(k,j,i)
541	local_pfi(2,i,j,k) = v(k,j,i)
542	local_pfi(3,i,j,k) = w(k,j,i)
543	tmp_norm = SQRT( u(k,j,i) * u(k,j,i) + &
544	v(k,j,i) * v(k,j,i) + &
545	w(k,j,i) * w(k,j,i) )
546	tmp_alpha = ACOS( 0.0_wp * u(k,j,i) / tmp_norm + &
547	0.0_wp * v(k,j,i) / tmp_norm - &
548	1.0_wp * w(k,j,i) / tmp_norm )
549	tmp_alpha_w = tmp_alpha / pi * 180.0_wp
550	local_pfi(4,i,j,k) = tmp_alpha_w
551	ENDDO
552	ENDDO
553	ENDDO
554
555	CALL cpu_log( log_point_s(31), 'dvrp_pathlines', 'start' )
556
557	CALL DVRP_DATA( m-1, local_pfi, 4, nx_dvrp, ny_dvrp, nz_dvrp, &
558	cyclic_dvrp, cyclic_dvrp, cyclic_dvrp )
559	CALL DVRP_VISUALIZE( m-1, 20, dvrp_filecount )
560
561	CALL cpu_log( log_point_s(31), 'dvrp_pathlines', 'stop' )
562
563	DEALLOCATE( local_pfi )
564
565	ENDIF
566
567	m = m + 1
568
569	ENDDO
570
571	dvrp_filecount = dvrp_filecount + 1
572
573	CALL cpu_log( log_point(27), 'data_output_dvrp', 'stop' )
574
575	#endif
576	END SUBROUTINE data_output_dvrp

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |