Home

Context Navigation

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

Last change on this file since 1310 was 1310, checked in by raasch, 10 years ago
update of GPL copyright
Property svn:keywords set to `Id`
File size: 17.6 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |