Home

Context Navigation

source: palm/trunk/SOURCE/data_output_mask.f90 @ 1034

Last change on this file since 1034 was 1031, checked in by raasch, 12 years ago
netCDF4 without parallel file support implemented additional define string netcdf4_parallel is required to switch on parallel file support
Property svn:keywords set to `Id`
File size: 16.7 KB

Rev	Line
[298]	1	SUBROUTINE data_output_mask( av )
	2
	3	!------------------------------------------------------------------------------!
[484]	4	! Current revisions:
[298]	5	! -----------------
[1031]	6	! netCDF4 without parallel file support implemented
[298]	7	!
	8	! Former revisions:
	9	! -----------------
	10	! $Id: data_output_mask.f90 1031 2012-10-19 14:35:30Z raasch $
	11	!
[1008]	12	! 1007 2012-09-19 14:30:36Z franke
	13	! Bugfix: calculation of pr must depend on the particle weighting factor,
	14	! missing calculation of ql_vp added
	15	!
[772]	16	! 771 2011-10-27 10:56:21Z heinze
	17	! +lpt
	18	!
[668]	19	! 667 2010-12-23 12:06:00Z suehring/gryschka
	20	! Calls of exchange_horiz are modified.
	21	!
[565]	22	! 564 2010-09-30 13:18:59Z helmke
	23	! start number of mask output files changed to 201, netcdf message identifiers
	24	! of masked output changed, palm message identifiers of masked output changed
	25	!
[494]	26	! 493 2010-03-01 08:30:24Z raasch
	27	! netcdf_format_mask* and format_parallel_io replaced by netcdf_data_format
	28	!
[482]	29	! 475 2010-02-04 02:26:16Z raasch
	30	! Bugfix in serial branch: arguments from array local_pf removed in N90_PUT_VAR
	31	!
[449]	32	! 410 2009-12-04 17:05:40Z letzel
	33	! Initial version
[298]	34	!
	35	! Description:
	36	! ------------
[1031]	37	! Masked data output in netCDF format for current mask (current value of mid).
[298]	38	!------------------------------------------------------------------------------!
	39
	40	#if defined( __netcdf )
	41	USE arrays_3d
	42	USE averaging
	43	USE cloud_parameters
	44	USE control_parameters
	45	USE cpulog
	46	USE grid_variables
	47	USE indices
	48	USE interfaces
	49	USE netcdf
	50	USE netcdf_control
	51	USE particle_attributes
	52	USE pegrid
	53
	54	IMPLICIT NONE
	55
	56	INTEGER :: av, ngp, file_id, i, if, is, j, k, l, n, psi, s, sender, &
	57	ind(6)
	58	LOGICAL :: found, resorted
	59	REAL :: mean_r, s_r3, s_r4
	60	REAL, DIMENSION(:,:,:), ALLOCATABLE :: local_pf
	61	#if defined( __parallel )
	62	REAL, DIMENSION(:,:,:), ALLOCATABLE :: total_pf
	63	#endif
	64	REAL, DIMENSION(:,:,:), POINTER :: to_be_resorted
	65
	66	!
	67	!-- Return, if nothing to output
	68	IF ( domask_no(mid,av) == 0 ) RETURN
	69
	70	CALL cpu_log (log_point(49),'data_output_mask','start')
	71
[493]	72	!
[298]	73	!-- Open output file.
[1031]	74	IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 4 ) ) &
[493]	75	THEN
[564]	76	CALL check_open( 200+mid+av*max_masks )
[409]	77	ENDIF
[298]	78
	79	!
	80	!-- Allocate total and local output arrays.
	81	#if defined( __parallel )
	82	IF ( myid == 0 ) THEN
	83	ALLOCATE( total_pf(mask_size(mid,1),mask_size(mid,2),mask_size(mid,3)) )
	84	ENDIF
	85	#endif
	86	ALLOCATE( local_pf(mask_size_l(mid,1),mask_size_l(mid,2), &
	87	mask_size_l(mid,3)) )
	88
	89	!
[1031]	90	!-- Update the netCDF time axis.
[298]	91	domask_time_count(mid,av) = domask_time_count(mid,av) + 1
[1031]	92	IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 4 ) ) &
[493]	93	THEN
[298]	94	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_time_mask(mid,av), &
	95	(/ simulated_time /), &
	96	start = (/ domask_time_count(mid,av) /), &
	97	count = (/ 1 /) )
[564]	98	CALL handle_netcdf_error( 'data_output_mask', 460 )
[298]	99	ENDIF
	100
	101	!
	102	!-- Loop over all variables to be written.
	103	if = 1
	104
	105	DO WHILE ( domask(mid,av,if)(1:1) /= ' ' )
	106	!
	107	!-- Reallocate local_pf on PE 0 since its shape changes during MPI exchange
[1031]	108	IF ( netcdf_data_format < 5 .AND. myid == 0 .AND. if > 1 ) THEN
[298]	109	DEALLOCATE( local_pf )
	110	ALLOCATE( local_pf(mask_size_l(mid,1),mask_size_l(mid,2), &
	111	mask_size_l(mid,3)) )
	112	ENDIF
	113	!
	114	!-- Store the variable chosen.
	115	resorted = .FALSE.
	116	SELECT CASE ( TRIM( domask(mid,av,if) ) )
	117
	118	CASE ( 'e' )
	119	IF ( av == 0 ) THEN
	120	to_be_resorted => e
	121	ELSE
	122	to_be_resorted => e_av
	123	ENDIF
	124
[771]	125	CASE ( 'lpt' )
	126	IF ( av == 0 ) THEN
	127	to_be_resorted => pt
	128	ELSE
	129	to_be_resorted => lpt_av
	130	ENDIF
	131
[298]	132	CASE ( 'p' )
	133	IF ( av == 0 ) THEN
	134	to_be_resorted => p
	135	ELSE
	136	to_be_resorted => p_av
	137	ENDIF
	138
	139	CASE ( 'pc' ) ! particle concentration (requires ghostpoint exchange)
	140	IF ( av == 0 ) THEN
	141	tend = prt_count
[667]	142	CALL exchange_horiz( tend, nbgp )
[298]	143	DO i = 1, mask_size_l(mid,1)
	144	DO j = 1, mask_size_l(mid,2)
	145	DO k = 1, mask_size_l(mid,3)
	146	local_pf(i,j,k) = tend(mask_k(mid,k), &
	147	mask_j(mid,j),mask_i(mid,i))
	148	ENDDO
	149	ENDDO
	150	ENDDO
	151	resorted = .TRUE.
	152	ELSE
[667]	153	CALL exchange_horiz( pc_av, nbgp )
[298]	154	to_be_resorted => pc_av
	155	ENDIF
	156
	157	CASE ( 'pr' ) ! mean particle radius
	158	IF ( av == 0 ) THEN
	159	DO i = nxl, nxr
	160	DO j = nys, nyn
	161	DO k = nzb, nzt+1
	162	psi = prt_start_index(k,j,i)
	163	s_r3 = 0.0
	164	s_r4 = 0.0
	165	DO n = psi, psi+prt_count(k,j,i)-1
[1007]	166	s_r3 = s_r3 + particles(n)%radius*3 &
	167	particles(n)%weight_factor
	168	s_r4 = s_r4 + particles(n)%radius*4 &
	169	particles(n)%weight_factor
[298]	170	ENDDO
	171	IF ( s_r3 /= 0.0 ) THEN
	172	mean_r = s_r4 / s_r3
	173	ELSE
	174	mean_r = 0.0
	175	ENDIF
	176	tend(k,j,i) = mean_r
	177	ENDDO
	178	ENDDO
	179	ENDDO
[667]	180	CALL exchange_horiz( tend, nbgp )
[298]	181	DO i = 1, mask_size_l(mid,1)
	182	DO j = 1, mask_size_l(mid,2)
	183	DO k = 1, mask_size_l(mid,3)
	184	local_pf(i,j,k) = tend(mask_k(mid,k), &
	185	mask_j(mid,j),mask_i(mid,i))
	186	ENDDO
	187	ENDDO
	188	ENDDO
	189	resorted = .TRUE.
	190	ELSE
[667]	191	CALL exchange_horiz( pr_av, nbgp )
[298]	192	to_be_resorted => pr_av
	193	ENDIF
	194
	195	CASE ( 'pt' )
	196	IF ( av == 0 ) THEN
	197	IF ( .NOT. cloud_physics ) THEN
	198	to_be_resorted => pt
	199	ELSE
	200	DO i = 1, mask_size_l(mid,1)
	201	DO j = 1, mask_size_l(mid,2)
	202	DO k = 1, mask_size_l(mid,3)
	203	local_pf(i,j,k) = &
	204	pt(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i)) &
	205	+ l_d_cp * pt_d_t(mask_k(mid,k)) * &
	206	ql(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i))
	207	ENDDO
	208	ENDDO
	209	ENDDO
	210	resorted = .TRUE.
	211	ENDIF
	212	ELSE
	213	to_be_resorted => pt_av
	214	ENDIF
	215
	216	CASE ( 'q' )
	217	IF ( av == 0 ) THEN
	218	to_be_resorted => q
	219	ELSE
	220	to_be_resorted => q_av
	221	ENDIF
[1007]	222
[298]	223	CASE ( 'ql' )
	224	IF ( av == 0 ) THEN
	225	to_be_resorted => ql
	226	ELSE
	227	to_be_resorted => ql_av
	228	ENDIF
	229
	230	CASE ( 'ql_c' )
	231	IF ( av == 0 ) THEN
	232	to_be_resorted => ql_c
	233	ELSE
	234	to_be_resorted => ql_c_av
	235	ENDIF
	236
	237	CASE ( 'ql_v' )
	238	IF ( av == 0 ) THEN
	239	to_be_resorted => ql_v
	240	ELSE
	241	to_be_resorted => ql_v_av
	242	ENDIF
	243
	244	CASE ( 'ql_vp' )
	245	IF ( av == 0 ) THEN
[1007]	246	DO i = nxl, nxr
	247	DO j = nys, nyn
	248	DO k = nzb, nz_do3d
	249	psi = prt_start_index(k,j,i)
	250	DO n = psi, psi+prt_count(k,j,i)-1
	251	tend(k,j,i) = tend(k,j,i) + &
	252	particles(n)%weight_factor / &
	253	prt_count(k,j,i)
	254	ENDDO
	255	ENDDO
	256	ENDDO
	257	ENDDO
	258	CALL exchange_horiz( tend, nbgp )
	259	DO i = 1, mask_size_l(mid,1)
	260	DO j = 1, mask_size_l(mid,2)
	261	DO k = 1, mask_size_l(mid,3)
	262	local_pf(i,j,k) = tend(mask_k(mid,k), &
	263	mask_j(mid,j),mask_i(mid,i))
	264	ENDDO
	265	ENDDO
	266	ENDDO
	267	resorted = .TRUE.
[298]	268	ELSE
[1007]	269	CALL exchange_horiz( ql_vp_av, nbgp )
[298]	270	to_be_resorted => ql_vp_av
	271	ENDIF
	272
	273	CASE ( 'qv' )
	274	IF ( av == 0 ) THEN
	275	DO i = 1, mask_size_l(mid,1)
	276	DO j = 1, mask_size_l(mid,2)
	277	DO k = 1, mask_size_l(mid,3)
	278	local_pf(i,j,k) = &
	279	q(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i)) - &
	280	ql(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i))
	281	ENDDO
	282	ENDDO
	283	ENDDO
	284	resorted = .TRUE.
	285	ELSE
	286	to_be_resorted => qv_av
	287	ENDIF
	288
	289	CASE ( 'rho' )
	290	IF ( av == 0 ) THEN
	291	to_be_resorted => rho
	292	ELSE
	293	to_be_resorted => rho_av
	294	ENDIF
[1007]	295
[298]	296	CASE ( 's' )
	297	IF ( av == 0 ) THEN
	298	to_be_resorted => q
	299	ELSE
[356]	300	to_be_resorted => s_av
[298]	301	ENDIF
[1007]	302
[298]	303	CASE ( 'sa' )
	304	IF ( av == 0 ) THEN
	305	to_be_resorted => sa
	306	ELSE
	307	to_be_resorted => sa_av
	308	ENDIF
[1007]	309
[298]	310	CASE ( 'u' )
	311	IF ( av == 0 ) THEN
	312	to_be_resorted => u
	313	ELSE
	314	to_be_resorted => u_av
	315	ENDIF
	316
	317	CASE ( 'v' )
	318	IF ( av == 0 ) THEN
	319	to_be_resorted => v
	320	ELSE
	321	to_be_resorted => v_av
	322	ENDIF
	323
	324	CASE ( 'vpt' )
	325	IF ( av == 0 ) THEN
	326	to_be_resorted => vpt
	327	ELSE
	328	to_be_resorted => vpt_av
	329	ENDIF
	330
	331	CASE ( 'w' )
	332	IF ( av == 0 ) THEN
	333	to_be_resorted => w
	334	ELSE
	335	to_be_resorted => w_av
	336	ENDIF
	337
	338	CASE DEFAULT
	339	!
	340	!-- User defined quantity
	341	CALL user_data_output_mask(av, domask(mid,av,if), found, local_pf )
	342	resorted = .TRUE.
	343
	344	IF ( .NOT. found ) THEN
	345	WRITE ( message_string, * ) 'no output available for: ', &
	346	TRIM( domask(mid,av,if) )
[564]	347	CALL message( 'data_output_mask', 'PA0327', 0, 0, 0, 6, 0 )
[298]	348	ENDIF
	349
	350	END SELECT
	351
	352	!
	353	!-- Resort the array to be output, if not done above
	354	IF ( .NOT. resorted ) THEN
	355	DO i = 1, mask_size_l(mid,1)
	356	DO j = 1, mask_size_l(mid,2)
	357	DO k = 1, mask_size_l(mid,3)
	358	local_pf(i,j,k) = to_be_resorted(mask_k(mid,k), &
	359	mask_j(mid,j),mask_i(mid,i))
	360	ENDDO
	361	ENDDO
	362	ENDDO
	363	ENDIF
	364
	365	!
	366	!-- I/O block. I/O methods are implemented
	367	!-- (1) for parallel execution
[1031]	368	!-- a. with netCDF 4 parallel I/O-enabled library
	369	!-- b. with netCDF 3 library
[298]	370	!-- (2) for serial execution.
	371	!-- The choice of method depends on the correct setting of preprocessor
[1031]	372	!-- directives __parallel and __netcdf4_parallel as well as on the parameter
[493]	373	!-- netcdf_data_format.
[298]	374	#if defined( __parallel )
[1031]	375	#if defined( __netcdf4_parallel )
	376	IF ( netcdf_data_format > 4 ) THEN
[298]	377	!
[1031]	378	!-- (1) a. Parallel I/O using netCDF 4 (not yet tested)
[298]	379	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), &
[340]	380	id_var_domask(mid,av,if), &
[409]	381	local_pf, &
[340]	382	start = (/ mask_start_l(mid,1), mask_start_l(mid,2), &
[409]	383	mask_start_l(mid,3), domask_time_count(mid,av) /), &
	384	count = (/ mask_size_l(mid,1), mask_size_l(mid,2), &
	385	mask_size_l(mid,3), 1 /) )
[564]	386	CALL handle_netcdf_error( 'data_output_mask', 461 )
[298]	387	ELSE
	388	#endif
	389	!
	390	!-- (1) b. Conventional I/O only through PE0
	391	!-- PE0 receives partial arrays from all processors of the respective mask
	392	!-- and outputs them. Here a barrier has to be set, because otherwise
	393	!-- "-MPI- FATAL: Remote protocol queue full" may occur.
	394	CALL MPI_BARRIER( comm2d, ierr )
	395
	396	ngp = mask_size_l(mid,1) * mask_size_l(mid,2) * mask_size_l(mid,3)
	397	IF ( myid == 0 ) THEN
	398	!
	399	!-- Local array can be relocated directly.
	400	total_pf( &
	401	mask_start_l(mid,1):mask_start_l(mid,1)+mask_size_l(mid,1)-1, &
	402	mask_start_l(mid,2):mask_start_l(mid,2)+mask_size_l(mid,2)-1, &
	403	mask_start_l(mid,3):mask_start_l(mid,3)+mask_size_l(mid,3)-1 ) &
	404	= local_pf
	405	!
	406	!-- Receive data from all other PEs.
	407	DO n = 1, numprocs-1
	408	!
	409	!-- Receive index limits first, then array.
	410	!-- Index limits are received in arbitrary order from the PEs.
	411	CALL MPI_RECV( ind(1), 6, MPI_INTEGER, MPI_ANY_SOURCE, 0, &
	412	comm2d, status, ierr )
	413	!
	414	!-- Not all PEs have data for the mask
	415	IF ( ind(1) /= -9999 ) THEN
	416	ngp = ( ind(2)-ind(1)+1 ) * (ind(4)-ind(3)+1 ) * &
	417	( ind(6)-ind(5)+1 )
	418	sender = status(MPI_SOURCE)
	419	DEALLOCATE( local_pf )
	420	ALLOCATE(local_pf(ind(1):ind(2),ind(3):ind(4),ind(5):ind(6)))
	421	CALL MPI_RECV( local_pf(ind(1),ind(3),ind(5)), ngp, &
	422	MPI_REAL, sender, 1, comm2d, status, ierr )
	423	total_pf(ind(1):ind(2),ind(3):ind(4),ind(5):ind(6)) &
	424	= local_pf
	425	ENDIF
	426	ENDDO
	427
	428	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), &
	429	id_var_domask(mid,av,if), total_pf, &
	430	start = (/ 1, 1, 1, domask_time_count(mid,av) /), &
	431	count = (/ mask_size(mid,1), mask_size(mid,2), &
	432	mask_size(mid,3), 1 /) )
[564]	433	CALL handle_netcdf_error( 'data_output_mask', 462 )
[298]	434
	435	ELSE
	436	!
	437	!-- If at least part of the mask resides on the PE, send the index
	438	!-- limits for the target array, otherwise send -9999 to PE0.
	439	IF ( mask_size_l(mid,1) > 0 .AND. mask_size_l(mid,2) > 0 .AND. &
	440	mask_size_l(mid,3) > 0 ) &
	441	THEN
	442	ind(1) = mask_start_l(mid,1)
	443	ind(2) = mask_start_l(mid,1) + mask_size_l(mid,1) - 1
	444	ind(3) = mask_start_l(mid,2)
	445	ind(4) = mask_start_l(mid,2) + mask_size_l(mid,2) - 1
	446	ind(5) = mask_start_l(mid,3)
	447	ind(6) = mask_start_l(mid,3) + mask_size_l(mid,3) - 1
	448	ELSE
	449	ind(1) = -9999; ind(2) = -9999
	450	ind(3) = -9999; ind(4) = -9999
	451	ind(5) = -9999; ind(6) = -9999
	452	ENDIF
	453	CALL MPI_SEND( ind(1), 6, MPI_INTEGER, 0, 0, comm2d, ierr )
	454	!
	455	!-- If applicable, send data to PE0.
	456	IF ( ind(1) /= -9999 ) THEN
	457	CALL MPI_SEND( local_pf(1,1,1), ngp, MPI_REAL, 0, 1, comm2d, &
	458	ierr )
	459	ENDIF
	460	ENDIF
	461	!
	462	!-- A barrier has to be set, because otherwise some PEs may proceed too
	463	!-- fast so that PE0 may receive wrong data on tag 0.
	464	CALL MPI_BARRIER( comm2d, ierr )
[1031]	465	#if defined( __netcdf4_parallel )
[298]	466	ENDIF
	467	#endif
	468	#else
	469	!
	470	!-- (2) For serial execution of PALM, the single processor (PE0) holds all
	471	!-- data and writes them directly to file.
	472	nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), &
	473	id_var_domask(mid,av,if), &
[475]	474	local_pf, &
[298]	475	start = (/ 1, 1, 1, domask_time_count(mid,av) /), &
	476	count = (/ mask_size_l(mid,1), mask_size_l(mid,2), &
	477	mask_size_l(mid,3), 1 /) )
[564]	478	CALL handle_netcdf_error( 'data_output_mask', 463 )
[298]	479	#endif
	480
	481	if = if + 1
[667]	482
[298]	483	ENDDO
	484
	485	!
	486	!-- Deallocate temporary arrays.
	487	DEALLOCATE( local_pf )
	488	#if defined( __parallel )
	489	IF ( myid == 0 ) THEN
	490	DEALLOCATE( total_pf )
	491	ENDIF
	492	#endif
	493
	494
	495	CALL cpu_log (log_point(49),'data_output_mask','stop','nobarrier')
	496	#endif
	497
	498	END SUBROUTINE data_output_mask

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |