Home

Context Navigation

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

Last change on this file since 1355 was 1354, checked in by heinze, 11 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 19.7 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |