Home

Context Navigation

source: palm/trunk/SOURCE/data_output_ptseries.f90 @ 1025

Last change on this file since 1025 was 826, checked in by raasch, 13 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 12.7 KB

Rev	Line
[1]	1	SUBROUTINE data_output_ptseries
	2
	3	!------------------------------------------------------------------------------!
[263]	4	! Current revisions:
[1]	5	! -----------------
[826]	6	!
[1]	7	!
	8	! Former revisions:
	9	! -----------------
[3]	10	! $Id: data_output_ptseries.f90 826 2012-02-19 03:41:34Z letzel $
[77]	11	!
[826]	12	! 825 2012-02-19 03:03:44Z raasch
	13	! mean/minimum/maximum particle radius added as output quantity,
	14	! particle attributes speed_x\|y\|z_sgs renamed rvar1\|2\|3
	15	!
[623]	16	! 622 2010-12-10 08:08:13Z raasch
	17	! optional barriers included in order to speed up collective operations
	18	!
[392]	19	! 291 2009-04-16 12:07:26Z raasch
	20	! simulated_time in NetCDF output replaced by time_since_reference_point.
	21	! Output of NetCDF messages with aid of message handling routine.
	22	!
[77]	23	! 60 2007-03-11 11:50:04Z raasch
	24	! Particles-package is now part of the default code.
	25	!
[3]	26	! RCS Log replace by Id keyword, revision history cleaned up
	27	!
[1]	28	! Revision 1.2 2006/08/22 13:51:13 raasch
	29	! Seperate output for particle groups
	30	!
	31	! Revision 1.1 2006/08/04 14:24:18 raasch
	32	! Initial revision
	33	!
	34	!
	35	! Description:
	36	! ------------
	37	! Output of particle data timeseries in NetCDF format.
	38	!------------------------------------------------------------------------------!
	39
[825]	40	USE cloud_parameters
[1]	41	USE control_parameters
	42	USE cpulog
	43	USE indices
	44	USE interfaces
	45	USE netcdf_control
	46	USE particle_attributes
	47	USE pegrid
	48
	49	IMPLICIT NONE
	50
	51
	52	INTEGER :: i, inum, j, n
	53
[825]	54	REAL, DIMENSION(:,:), ALLOCATABLE :: pts_value, pts_value_l
[1]	55
	56
	57	CALL cpu_log( log_point(36), 'data_output_ptseries', 'start' )
	58
	59	IF ( myid == 0 .AND. netcdf_output ) THEN
	60	!
	61	!-- Open file for time series output in NetCDF format
	62	dopts_time_count = dopts_time_count + 1
	63	CALL check_open( 109 )
	64	#if defined( __netcdf )
	65	!
	66	!-- Update the particle time series time axis
[291]	67	nc_stat = NF90_PUT_VAR( id_set_pts, id_var_time_pts, &
	68	(/ time_since_reference_point /), &
[1]	69	start = (/ dopts_time_count /), count = (/ 1 /) )
[263]	70	CALL handle_netcdf_error( 'data_output_ptseries', 391 )
[1]	71	#endif
	72
	73	ENDIF
	74
[825]	75	ALLOCATE( pts_value(0:number_of_particle_groups,dopts_num), &
	76	pts_value_l(0:number_of_particle_groups,dopts_num) )
	77
[1]	78	pts_value_l = 0.0
[825]	79	pts_value_l(:,16) = 9999999.9 ! for calculation of minimum radius
[1]	80
	81	!
	82	!-- Calculate or collect the particle time series quantities for all particles
	83	!-- and seperately for each particle group (if there is more than one group)
	84	DO n = 1, number_of_particles
	85
	86	pts_value_l(0,1) = number_of_particles ! total # of particles
	87	pts_value_l(0,2) = pts_value_l(0,2) + &
	88	( particles(n)%x - particles(n)%origin_x ) ! mean x
	89	pts_value_l(0,3) = pts_value_l(0,3) + &
	90	( particles(n)%y - particles(n)%origin_y ) ! mean y
	91	pts_value_l(0,4) = pts_value_l(0,4) + &
	92	( particles(n)%z - particles(n)%origin_z ) ! mean z
	93	pts_value_l(0,5) = pts_value_l(0,5) + particles(n)%z ! mean z (absolute)
	94	pts_value_l(0,6) = pts_value_l(0,6) + particles(n)%speed_x ! mean u
	95	pts_value_l(0,7) = pts_value_l(0,7) + particles(n)%speed_y ! mean v
	96	pts_value_l(0,8) = pts_value_l(0,8) + particles(n)%speed_z ! mean w
[825]	97	IF ( .NOT. curvature_solution_effects ) THEN
	98	pts_value_l(0,9) = pts_value_l(0,9) + particles(n)%rvar1 ! mean sgsu
	99	pts_value_l(0,10) = pts_value_l(0,10) + particles(n)%rvar2 ! mean sgsv
	100	pts_value_l(0,11) = pts_value_l(0,11) + particles(n)%rvar3 ! mean sgsw
	101	ENDIF
[1]	102	IF ( particles(n)%speed_z > 0.0 ) THEN
	103	pts_value_l(0,12) = pts_value_l(0,12) + 1.0 ! # of upward moving prts
	104	pts_value_l(0,13) = pts_value_l(0,13) + &
	105	particles(n)%speed_z ! mean w upw.
	106	ELSE
	107	pts_value_l(0,14) = pts_value_l(0,14) + &
	108	particles(n)%speed_z ! mean w down
	109	ENDIF
[825]	110	pts_value_l(0,15) = pts_value_l(0,15) + particles(n)%radius ! mean rad
	111	pts_value_l(0,16) = MIN( pts_value_l(0,16), particles(n)%radius ) ! minrad
	112	pts_value_l(0,17) = MAX( pts_value_l(0,17), particles(n)%radius ) ! maxrad
	113	pts_value_l(0,18) = number_of_particles
	114	pts_value_l(0,19) = number_of_particles
[1]	115
	116	!
	117	!-- Repeat the same for the respective particle group
	118	IF ( number_of_particle_groups > 1 ) THEN
	119	j = particles(n)%group
	120
	121	pts_value_l(j,1) = pts_value_l(j,1) + 1
	122	pts_value_l(j,2) = pts_value_l(j,2) + &
	123	( particles(n)%x - particles(n)%origin_x )
	124	pts_value_l(j,3) = pts_value_l(j,3) + &
	125	( particles(n)%y - particles(n)%origin_y )
	126	pts_value_l(j,4) = pts_value_l(j,4) + &
	127	( particles(n)%z - particles(n)%origin_z )
	128	pts_value_l(j,5) = pts_value_l(j,5) + particles(n)%z
	129	pts_value_l(j,6) = pts_value_l(j,6) + particles(n)%speed_x
	130	pts_value_l(j,7) = pts_value_l(j,7) + particles(n)%speed_y
	131	pts_value_l(j,8) = pts_value_l(j,8) + particles(n)%speed_z
[825]	132	IF ( .NOT. curvature_solution_effects ) THEN
	133	pts_value_l(j,9) = pts_value_l(j,9) + particles(n)%rvar1
	134	pts_value_l(j,10) = pts_value_l(j,10) + particles(n)%rvar2
	135	pts_value_l(j,11) = pts_value_l(j,11) + particles(n)%rvar3
	136	ENDIF
[1]	137	IF ( particles(n)%speed_z > 0.0 ) THEN
	138	pts_value_l(j,12) = pts_value_l(j,12) + 1.0
	139	pts_value_l(j,13) = pts_value_l(j,13) + particles(n)%speed_z
	140	ELSE
	141	pts_value_l(j,14) = pts_value_l(j,14) + particles(n)%speed_z
	142	ENDIF
[824]	143	pts_value_l(j,15) = pts_value_l(j,15) + particles(n)%radius
[825]	144	pts_value_l(j,16) = MIN( pts_value(j,16), particles(n)%radius )
	145	pts_value_l(j,17) = MAX( pts_value(j,17), particles(n)%radius )
	146	pts_value_l(j,18) = pts_value_l(j,18) + 1.0
	147	pts_value_l(j,19) = pts_value_l(j,19) + 1.0
[1]	148
	149	ENDIF
	150
	151	ENDDO
	152
	153	#if defined( __parallel )
	154	!
	155	!-- Sum values of the subdomains
	156	inum = number_of_particle_groups + 1
	157
[622]	158	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[824]	159	CALL MPI_ALLREDUCE( pts_value_l(0,1), pts_value(0,1), 15*inum, MPI_REAL, &
[1]	160	MPI_SUM, comm2d, ierr )
[622]	161	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[824]	162	CALL MPI_ALLREDUCE( pts_value_l(0,16), pts_value(0,16), inum, MPI_REAL, &
[825]	163	MPI_MIN, comm2d, ierr )
	164	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	165	CALL MPI_ALLREDUCE( pts_value_l(0,17), pts_value(0,17), inum, MPI_REAL, &
[1]	166	MPI_MAX, comm2d, ierr )
[622]	167	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[825]	168	CALL MPI_ALLREDUCE( pts_value_l(0,18), pts_value(0,18), inum, MPI_REAL, &
	169	MPI_MAX, comm2d, ierr )
	170	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	171	CALL MPI_ALLREDUCE( pts_value_l(0,19), pts_value(0,19), inum, MPI_REAL, &
[1]	172	MPI_MIN, comm2d, ierr )
	173	#else
[825]	174	pts_value(:,1:19) = pts_value_l(:,1:19)
[1]	175	#endif
	176
	177	!
[825]	178	!-- Normalize the above calculated quantities (except min/max values) with the
	179	!-- total number of particles
[1]	180	IF ( number_of_particle_groups > 1 ) THEN
	181	inum = number_of_particle_groups
	182	ELSE
	183	inum = 0
	184	ENDIF
	185
	186	DO j = 0, inum
	187
	188	IF ( pts_value(j,1) > 0.0 ) THEN
	189
[824]	190	pts_value(j,2:15) = pts_value(j,2:15) / pts_value(j,1)
[1]	191	IF ( pts_value(j,12) > 0.0 .AND. pts_value(j,12) < 1.0 ) THEN
	192	pts_value(j,13) = pts_value(j,13) / pts_value(j,12)
	193	pts_value(j,14) = pts_value(j,14) / ( 1.0 - pts_value(j,12) )
	194	ELSEIF ( pts_value(j,12) == 0.0 ) THEN
	195	pts_value(j,13) = -1.0
	196	ELSE
	197	pts_value(j,14) = -1.0
	198	ENDIF
	199
	200	ENDIF
	201
	202	ENDDO
	203
	204	!
	205	!-- Calculate higher order moments of particle time series quantities,
	206	!-- seperately for each particle group (if there is more than one group)
	207	DO n = 1, number_of_particles
	208
[825]	209	pts_value_l(0,20) = pts_value_l(0,20) + ( particles(n)%x - &
[1]	210	particles(n)%origin_x - pts_value(0,2) )*2 ! x2
[825]	211	pts_value_l(0,21) = pts_value_l(0,21) + ( particles(n)%y - &
[1]	212	particles(n)%origin_y - pts_value(0,3) )*2 ! y2
[825]	213	pts_value_l(0,22) = pts_value_l(0,22) + ( particles(n)%z - &
[1]	214	particles(n)%origin_z - pts_value(0,4) )*2 ! z2
[825]	215	pts_value_l(0,23) = pts_value_l(0,23) + ( particles(n)%speed_x - &
	216	pts_value(0,6) )*2 ! u2
	217	pts_value_l(0,24) = pts_value_l(0,24) + ( particles(n)%speed_y - &
	218	pts_value(0,7) )*2 ! v2
	219	pts_value_l(0,25) = pts_value_l(0,25) + ( particles(n)%speed_z - &
	220	pts_value(0,8) )*2 ! w2
	221	IF ( .NOT. curvature_solution_effects ) THEN
	222	pts_value_l(0,26) = pts_value_l(0,26) + ( particles(n)%rvar1 - &
	223	pts_value(0,9) )**2 ! u"2
	224	pts_value_l(0,27) = pts_value_l(0,27) + ( particles(n)%rvar2 - &
	225	pts_value(0,10) )**2 ! v"2
	226	pts_value_l(0,28) = pts_value_l(0,28) + ( particles(n)%rvar3 - &
	227	pts_value(0,11) )**2 ! w"2
	228	ENDIF
[1]	229	!
	230	!-- Repeat the same for the respective particle group
	231	IF ( number_of_particle_groups > 1 ) THEN
	232	j = particles(n)%group
	233
[825]	234	pts_value_l(j,20) = pts_value_l(j,20) + ( particles(n)%x - &
[1]	235	particles(n)%origin_x - pts_value(j,2) )**2
[825]	236	pts_value_l(j,21) = pts_value_l(j,21) + ( particles(n)%y - &
[1]	237	particles(n)%origin_y - pts_value(j,3) )**2
[825]	238	pts_value_l(j,22) = pts_value_l(j,22) + ( particles(n)%z - &
[1]	239	particles(n)%origin_z - pts_value(j,4) )**2
[825]	240	pts_value_l(j,23) = pts_value_l(j,23) + ( particles(n)%speed_x - &
	241	pts_value(j,6) )**2
	242	pts_value_l(j,24) = pts_value_l(j,24) + ( particles(n)%speed_y - &
	243	pts_value(j,7) )**2
	244	pts_value_l(j,25) = pts_value_l(j,25) + ( particles(n)%speed_z - &
	245	pts_value(j,8) )**2
	246	IF ( .NOT. curvature_solution_effects ) THEN
	247	pts_value_l(j,26) = pts_value_l(j,26) + ( particles(n)%rvar1 - &
	248	pts_value(j,9) )**2
	249	pts_value_l(j,27) = pts_value_l(j,27) + ( particles(n)%rvar2 - &
	250	pts_value(j,10) )**2
	251	pts_value_l(j,28) = pts_value_l(j,28) + ( particles(n)%rvar3 - &
	252	pts_value(j,11) )**2
	253	ENDIF
[1]	254	ENDIF
	255
	256	ENDDO
	257
[825]	258	pts_value_l(0,29) = ( number_of_particles - pts_value(0,1) / numprocs )**2
[1]	259	! variance of particle numbers
	260	IF ( number_of_particle_groups > 1 ) THEN
	261	DO j = 1, number_of_particle_groups
[825]	262	pts_value_l(j,29) = ( pts_value_l(j,1) - &
[1]	263	pts_value(j,1) / numprocs )**2
	264	ENDDO
	265	ENDIF
	266
	267	#if defined( __parallel )
	268	!
	269	!-- Sum values of the subdomains
	270	inum = number_of_particle_groups + 1
	271
[622]	272	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
[825]	273	CALL MPI_ALLREDUCE( pts_value_l(0,20), pts_value(0,20), inum*10, MPI_REAL, &
[1]	274	MPI_SUM, comm2d, ierr )
	275	#else
[825]	276	pts_value(:,20:29) = pts_value_l(:,20:29)
[1]	277	#endif
	278
	279	!
	280	!-- Normalize the above calculated quantities with the total number of
	281	!-- particles
	282	IF ( number_of_particle_groups > 1 ) THEN
	283	inum = number_of_particle_groups
	284	ELSE
	285	inum = 0
	286	ENDIF
	287
	288	DO j = 0, inum
	289
	290	IF ( pts_value(j,1) > 0.0 ) THEN
[825]	291	pts_value(j,20:28) = pts_value(j,20:28) / pts_value(j,1)
[1]	292	ENDIF
[825]	293	pts_value(j,29) = pts_value(j,29) / numprocs
[1]	294
	295	ENDDO
	296
	297	#if defined( __netcdf )
	298	!
	299	!-- Output particle time series quantities in NetCDF format
	300	IF ( myid == 0 .AND. netcdf_output ) THEN
	301	DO j = 0, inum
	302	DO i = 1, dopts_num
	303	nc_stat = NF90_PUT_VAR( id_set_pts, id_var_dopts(i,j), &
	304	(/ pts_value(j,i) /), &
	305	start = (/ dopts_time_count /), &
	306	count = (/ 1 /) )
[263]	307	CALL handle_netcdf_error( 'data_output_ptseries', 392 )
[1]	308	ENDDO
	309	ENDDO
	310	ENDIF
	311	#endif
	312
[825]	313	DEALLOCATE( pts_value, pts_value_l )
	314
[1]	315	CALL cpu_log( log_point(36), 'data_output_ptseries','stop', 'nobarrier' )
	316
	317	END SUBROUTINE data_output_ptseries

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |