Context Navigation

concoutput_mpi.f90 @ f3054ea

GFS_025dev

Last change on this file since f3054ea was f3054ea, checked in by Espen Sollum <eso@…>, 4 years ago
Changed from grib_api to eccodes. MPI: implemented linversionout=1; fixed calculation of grid_initial fields.
Property mode set to `100644`
File size: 26.7 KB

Rev	Line
[8a65cb0]	1	!**********************************************************************
	2	! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010 *
	3	! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa, *
	4	! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann *
	5	! *
	6	! This file is part of FLEXPART. *
	7	! *
	8	! FLEXPART is free software: you can redistribute it and/or modify *
	9	! it under the terms of the GNU General Public License as published by*
	10	! the Free Software Foundation, either version 3 of the License, or *
	11	! (at your option) any later version. *
	12	! *
	13	! FLEXPART is distributed in the hope that it will be useful, *
	14	! but WITHOUT ANY WARRANTY; without even the implied warranty of *
	15	! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
	16	! GNU General Public License for more details. *
	17	! *
	18	! You should have received a copy of the GNU General Public License *
	19	! along with FLEXPART. If not, see <http://www.gnu.org/licenses/>. *
	20	!**********************************************************************
	21
	22	subroutine concoutput(itime,outnum,gridtotalunc,wetgridtotalunc, &
	23	drygridtotalunc)
	24	! i i o o
	25	! o
	26	!*****************************************************************************
	27	! *
	28	! Output of the concentration grid and the receptor concentrations. *
	29	! *
	30	! Author: A. Stohl *
	31	! *
	32	! 24 May 1995 *
	33	! *
	34	! 13 April 1999, Major update: if output size is smaller, dump output *
	35	! in sparse matrix format; additional output of *
	36	! uncertainty *
	37	! *
	38	! 05 April 2000, Major update: output of age classes; output for backward*
	39	! runs is time spent in grid cell times total mass of *
	40	! species. *
	41	! *
	42	! 17 February 2002, Appropriate dimensions for backward and forward runs *
	43	! are now specified in file par_mod *
	44	! *
	45	! June 2006, write grid in sparse matrix with a single write command *
	46	! in order to save disk space *
	47	! *
	48	! 2008 new sparse matrix format *
	49	! *
	50	! Changes eso: *
	51	! 2014 MPI version This routine is only called by root MPI *
	52	! process (the other processes have sent *
	53	! their fields to root) *
	54	!*****************************************************************************
	55	! *
	56	! Variables: *
	57	! outnum number of samples *
	58	! ncells number of cells with non-zero concentrations *
	59	! sparse .true. if in sparse matrix format, else .false. *
	60	! tot_mu 1 for forward, initial mass mixing ration for backw. runs *
	61	! *
	62	!*****************************************************************************
	63
	64	use unc_mod
	65	use point_mod
	66	use outg_mod
	67	use par_mod
	68	use com_mod
	69	use mpi_mod
[6a678e3]	70	use mean_mod
[8a65cb0]	71
	72	implicit none
	73
	74	real(kind=dp) :: jul
	75	integer :: itime,i,ix,jy,kz,ks,kp,l,iix,jjy,kzz,nage,jjjjmmdd,ihmmss
	76	integer :: sp_count_i,sp_count_r
	77	real :: sp_fact
	78	real :: outnum,densityoutrecept(maxreceptor),xl,yl
[f3054ea]	79	! RLT
	80	real :: densitydryrecept(maxreceptor)
	81	real :: factor_dryrecept(maxreceptor)
[8a65cb0]	82
	83	!real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
	84	! +grid(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,maxpointspec_act,
	85	! + maxageclass)
	86	!real wetgrid(0:numxgrid-1,0:numygrid-1,maxspec,maxpointspec_act,
	87	! + maxageclass)
	88	!real drygrid(0:numxgrid-1,0:numygrid-1,maxspec,
	89	! + maxpointspec_act,maxageclass)
	90	!real gridsigma(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,
	91	! + maxpointspec_act,maxageclass),
	92	! + drygridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
	93	! + maxpointspec_act,maxageclass),
	94	! + wetgridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
	95	! + maxpointspec_act,maxageclass)
	96	!real factor(0:numxgrid-1,0:numygrid-1,numzgrid)
	97	!real sparse_dump_r(numxgridnumygridnumzgrid)
	98	!integer sparse_dump_i(numxgridnumygridnumzgrid)
	99
	100	!real sparse_dump_u(numxgridnumygridnumzgrid)
[6a678e3]	101	real(dep_prec) :: auxgrid(nclassunc)
	102	real(sp) :: gridtotal,gridsigmatotal,gridtotalunc
	103	real(dep_prec) :: wetgridtotal,wetgridsigmatotal,wetgridtotalunc
	104	real(dep_prec) :: drygridtotal,drygridsigmatotal,drygridtotalunc
[8a65cb0]	105	real :: halfheight,dz,dz1,dz2,tot_mu(maxspec,maxpointspec_act)
	106	real,parameter :: smallnum = tiny(0.0) ! smallest number that can be handled
	107	real,parameter :: weightair=28.97
	108	logical :: sp_zer
[20963b1]	109	logical,save :: init=.true.
[8a65cb0]	110	character :: adate8,atime6
	111	character(len=3) :: anspec
[5f9d14a]	112	integer :: mind
[b069789]	113	! mind eso: added to ensure identical results between 2&3-fields versions
	114	character(LEN=8),save :: file_stat='REPLACE'
	115	logical :: ldates_file
[f3054ea]	116	logical :: lexist
[b069789]	117	integer :: ierr
	118	character(LEN=100) :: dates_char
[8a65cb0]	119
	120	! Measure execution time
	121	if (mp_measure_time) call mpif_mtime('rootonly',0)
	122
	123	! Determine current calendar date, needed for the file name
	124	!**********************************************************
	125
	126	jul=bdate+real(itime,kind=dp)/86400._dp
	127	call caldate(jul,jjjjmmdd,ihmmss)
	128	write(adate,'(i8.8)') jjjjmmdd
	129	write(atime,'(i6.6)') ihmmss
[b069789]	130
	131	! Overwrite existing dates file on first call, later append to it
	132	! This fixes a bug where the dates file kept growing across multiple runs
	133
[20963b1]	134	! If 'dates' file exists in output directory, make a backup
[b069789]	135	inquire(file=path(2)(1:length(2))//'dates', exist=ldates_file)
	136	if (ldates_file.and.init) then
	137	open(unit=unitdates, file=path(2)(1:length(2))//'dates',form='formatted', &
	138	&access='sequential', status='old', action='read', iostat=ierr)
	139	open(unit=unittmp, file=path(2)(1:length(2))//'dates.bak', access='sequential', &
	140	&status='replace', action='write', form='formatted', iostat=ierr)
	141	do while (.true.)
	142	read(unitdates, '(a)', iostat=ierr) dates_char
	143	if (ierr.ne.0) exit
	144	write(unit=unittmp, fmt='(a)', iostat=ierr, advance='yes') trim(dates_char)
	145	! if (ierr.ne.0) write(,) "Write error, ", ierr
	146	end do
	147	close(unit=unitdates)
	148	close(unit=unittmp)
	149	end if
	150
	151	open(unitdates,file=path(2)(1:length(2))//'dates', ACCESS='APPEND', STATUS=file_stat)
[8a65cb0]	152	write(unitdates,'(a)') adate//atime
	153	close(unitdates)
	154
[478e9e6]	155	! Overwrite existing dates file on first call, later append to it
	156	! This fixes a bug where the dates file kept growing across multiple runs
	157	IF (init) THEN
	158	file_stat='OLD'
	159	init=.false.
	160	END IF
	161
	162
[8a65cb0]	163	! For forward simulations, output fields have dimension MAXSPEC,
	164	! for backward simulations, output fields have dimension MAXPOINT.
	165	! Thus, make loops either about nspec, or about numpoint
	166	!*****************************************************************
	167
	168
	169	if (ldirect.eq.1) then
	170	do ks=1,nspec
	171	do kp=1,maxpointspec_act
	172	tot_mu(ks,kp)=1
	173	end do
	174	end do
	175	else
	176	do ks=1,nspec
	177	do kp=1,maxpointspec_act
	178	tot_mu(ks,kp)=xmass(kp,ks)
	179	end do
	180	end do
	181	endif
	182
	183
	184	!*******************************************************************
	185	! Compute air density: sufficiently accurate to take it
	186	! from coarse grid at some time
	187	! Determine center altitude of output layer, and interpolate density
	188	! data to that altitude
	189	!*******************************************************************
	190
[5f9d14a]	191	mind=memind(2)
[8a65cb0]	192	do kz=1,numzgrid
	193	if (kz.eq.1) then
	194	halfheight=outheight(1)/2.
	195	else
	196	halfheight=(outheight(kz)+outheight(kz-1))/2.
	197	endif
	198	do kzz=2,nz
	199	if ((height(kzz-1).lt.halfheight).and. &
	200	(height(kzz).gt.halfheight)) goto 46
	201	end do
	202	46 kzz=max(min(kzz,nz),2)
	203	dz1=halfheight-height(kzz-1)
	204	dz2=height(kzz)-halfheight
	205	dz=dz1+dz2
	206	do jy=0,numygrid-1
	207	do ix=0,numxgrid-1
	208	xl=outlon0+real(ix)*dxout
	209	yl=outlat0+real(jy)*dyout
	210	xl=(xl-xlon0)/dx
	211	yl=(yl-ylat0)/dy !v9.1.1
	212	iix=max(min(nint(xl),nxmin1),0)
	213	jjy=max(min(nint(yl),nymin1),0)
[5f9d14a]	214	! densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
	215	! rho(iix,jjy,kzz-1,2)*dz2)/dz
	216	densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,mind)*dz1+ &
	217	rho(iix,jjy,kzz-1,mind)*dz2)/dz
[f3054ea]	218	! RLT
	219	densitydrygrid(ix,jy,kz)=(rho_dry(iix,jjy,kzz,mind)*dz1+ &
	220	rho_dry(iix,jjy,kzz-1,mind)*dz2)/dz
[8a65cb0]	221	end do
	222	end do
	223	end do
	224
	225	do i=1,numreceptor
	226	xl=xreceptor(i)
	227	yl=yreceptor(i)
	228	iix=max(min(nint(xl),nxmin1),0)
	229	jjy=max(min(nint(yl),nymin1),0)
[5f9d14a]	230	!densityoutrecept(i)=rho(iix,jjy,1,2)
	231	densityoutrecept(i)=rho(iix,jjy,1,mind)
[f3054ea]	232	! RLT
	233	densitydryrecept(i)=rho_dry(iix,jjy,1,mind)
[8a65cb0]	234	end do
	235
[f3054ea]	236	! RLT
	237	! conversion factor for output relative to dry air
	238	factor_drygrid=densityoutgrid/densitydrygrid
	239	factor_dryrecept=densityoutrecept/densitydryrecept
[8a65cb0]	240
	241	! Output is different for forward and backward simulations
	242	do kz=1,numzgrid
	243	do jy=0,numygrid-1
	244	do ix=0,numxgrid-1
	245	if (ldirect.eq.1) then
	246	factor3d(ix,jy,kz)=1.e12/volume(ix,jy,kz)/outnum
	247	else
	248	factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
	249	endif
	250	end do
	251	end do
	252	end do
	253
	254	!*********************************************************************
	255	! Determine the standard deviation of the mean concentration or mixing
	256	! ratio (uncertainty of the output) and the dry and wet deposition
	257	!*********************************************************************
	258
	259	gridtotal=0.
	260	gridsigmatotal=0.
	261	gridtotalunc=0.
	262	wetgridtotal=0.
	263	wetgridsigmatotal=0.
	264	wetgridtotalunc=0.
	265	drygridtotal=0.
	266	drygridsigmatotal=0.
	267	drygridtotalunc=0.
	268
	269	do ks=1,nspec
	270
	271	write(anspec,'(i3.3)') ks
	272
[20963b1]	273	if (DRYBKDEP.or.WETBKDEP) then !scavdep output
	274	if (DRYBKDEP) &
	275	open(unitoutgrid,file=path(2)(1:length(2))//'grid_drydep_'//adate// &
	276	atime//'_'//anspec,form='unformatted')
	277	if (WETBKDEP) &
	278	open(unitoutgrid,file=path(2)(1:length(2))//'grid_wetdep_'//adate// &
[8a65cb0]	279	atime//'_'//anspec,form='unformatted')
[20963b1]	280	write(unitoutgrid) itime
	281	else
	282	if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
	283	if (ldirect.eq.1) then
	284	open(unitoutgrid,file=path(2)(1:length(2))//'grid_conc_'//adate// &
	285	atime//'_'//anspec,form='unformatted')
	286	else
	287	open(unitoutgrid,file=path(2)(1:length(2))//'grid_time_'//adate// &
	288	atime//'_'//anspec,form='unformatted')
	289	endif
	290	write(unitoutgrid) itime
	291	endif
	292	if ((iout.eq.2).or.(iout.eq.3)) then ! mixing ratio
	293	open(unitoutgridppt,file=path(2)(1:length(2))//'grid_pptv_'//adate// &
	294	atime//'_'//anspec,form='unformatted')
	295	write(unitoutgridppt) itime
	296	endif
	297	endif ! if deposition output
[8a65cb0]	298
	299	do kp=1,maxpointspec_act
	300	do nage=1,nageclass
	301
	302	do jy=0,numygrid-1
	303	do ix=0,numxgrid-1
	304
	305	! WET DEPOSITION
	306	if ((WETDEP).and.(ldirect.gt.0)) then
	307	do l=1,nclassunc
	308	auxgrid(l)=wetgridunc0(ix,jy,ks,kp,l,nage)
	309	end do
	310	call mean(auxgrid,wetgrid(ix,jy), &
	311	wetgridsigma(ix,jy),nclassunc)
	312	! Multiply by number of classes to get total concentration
	313	wetgrid(ix,jy)=wetgrid(ix,jy) &
	314	*nclassunc
	315	wetgridtotal=wetgridtotal+wetgrid(ix,jy)
	316	! Calculate standard deviation of the mean
	317	wetgridsigma(ix,jy)= &
	318	wetgridsigma(ix,jy)* &
	319	sqrt(real(nclassunc))
	320	wetgridsigmatotal=wetgridsigmatotal+ &
	321	wetgridsigma(ix,jy)
	322	endif
	323
	324	! DRY DEPOSITION
	325	if ((DRYDEP).and.(ldirect.gt.0)) then
	326	do l=1,nclassunc
	327	auxgrid(l)=drygridunc0(ix,jy,ks,kp,l,nage)
	328	end do
	329	call mean(auxgrid,drygrid(ix,jy), &
	330	drygridsigma(ix,jy),nclassunc)
	331	! Multiply by number of classes to get total concentration
	332	drygrid(ix,jy)=drygrid(ix,jy)* &
	333	nclassunc
	334	drygridtotal=drygridtotal+drygrid(ix,jy)
	335	! Calculate standard deviation of the mean
	336	drygridsigma(ix,jy)= &
	337	drygridsigma(ix,jy)* &
	338	sqrt(real(nclassunc))
	339	drygridsigmatotal=drygridsigmatotal+ &
	340	drygridsigma(ix,jy)
	341	endif
	342
	343	! CONCENTRATION OR MIXING RATIO
	344	do kz=1,numzgrid
	345	do l=1,nclassunc
	346	auxgrid(l)=gridunc(ix,jy,kz,ks,kp,l,nage)
	347	end do
	348	call mean(auxgrid,grid(ix,jy,kz), &
	349	gridsigma(ix,jy,kz),nclassunc)
	350	! Multiply by number of classes to get total concentration
	351	grid(ix,jy,kz)= &
	352	grid(ix,jy,kz)*nclassunc
	353	gridtotal=gridtotal+grid(ix,jy,kz)
	354	! Calculate standard deviation of the mean
	355	gridsigma(ix,jy,kz)= &
	356	gridsigma(ix,jy,kz)* &
	357	sqrt(real(nclassunc))
	358	gridsigmatotal=gridsigmatotal+ &
	359	gridsigma(ix,jy,kz)
	360	end do
	361	end do
	362	end do
	363
	364
	365
	366
	367	!*******************************************************************
	368	! Generate output: may be in concentration (ng/m3) or in mixing
	369	! ratio (ppt) or both
	370	! Output the position and the values alternated multiplied by
	371	! 1 or -1, first line is number of values, number of positions
	372	! For backward simulations, the unit is seconds, stored in grid_time
	373	!*******************************************************************
	374
	375	! Concentration output
	376	!*********************
[6741557]	377	if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
[8a65cb0]	378
	379	! Wet deposition
	380	sp_count_i=0
	381	sp_count_r=0
	382	sp_fact=-1.
	383	sp_zer=.true.
	384	if ((ldirect.eq.1).and.(WETDEP)) then
	385	do jy=0,numygrid-1
	386	do ix=0,numxgrid-1
	387	!oncentraion greater zero
	388	if (wetgrid(ix,jy).gt.smallnum) then
	389	if (sp_zer.eqv..true.) then ! first non zero value
	390	sp_count_i=sp_count_i+1
	391	sparse_dump_i(sp_count_i)=ix+jy*numxgrid
	392	sp_zer=.false.
	393	sp_fact=sp_fact*(-1.)
	394	endif
	395	sp_count_r=sp_count_r+1
	396	sparse_dump_r(sp_count_r)= &
	397	sp_fact1.e12wetgrid(ix,jy)/area(ix,jy)
	398	! sparse_dump_u(sp_count_r)=
	399	!+ 1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
	400	else ! concentration is zero
	401	sp_zer=.true.
	402	endif
	403	end do
	404	end do
	405	else
	406	sp_count_i=0
	407	sp_count_r=0
	408	endif
	409	write(unitoutgrid) sp_count_i
	410	write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
	411	write(unitoutgrid) sp_count_r
	412	write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
	413	! write(unitoutgrid) sp_count_u
	414	! write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
	415
	416	! Dry deposition
	417	sp_count_i=0
	418	sp_count_r=0
	419	sp_fact=-1.
	420	sp_zer=.true.
	421	if ((ldirect.eq.1).and.(DRYDEP)) then
	422	do jy=0,numygrid-1
	423	do ix=0,numxgrid-1
	424	if (drygrid(ix,jy).gt.smallnum) then
	425	if (sp_zer.eqv..true.) then ! first non zero value
	426	sp_count_i=sp_count_i+1
	427	sparse_dump_i(sp_count_i)=ix+jy*numxgrid
	428	sp_zer=.false.
	429	sp_fact=sp_fact*(-1.)
	430	endif
	431	sp_count_r=sp_count_r+1
	432	sparse_dump_r(sp_count_r)= &
	433	sp_fact* &
	434	1.e12*drygrid(ix,jy)/area(ix,jy)
	435	! sparse_dump_u(sp_count_r)=
	436	!+ 1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
	437	else ! concentration is zero
	438	sp_zer=.true.
	439	endif
	440	end do
	441	end do
	442	else
	443	sp_count_i=0
	444	sp_count_r=0
	445	endif
	446	write(unitoutgrid) sp_count_i
	447	write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
	448	write(unitoutgrid) sp_count_r
	449	write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
	450	! write(,) sp_count_u
	451	! write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
	452
	453
	454
	455	! Concentrations
	456	sp_count_i=0
	457	sp_count_r=0
	458	sp_fact=-1.
	459	sp_zer=.true.
	460	do kz=1,numzgrid
	461	do jy=0,numygrid-1
	462	do ix=0,numxgrid-1
	463	if (grid(ix,jy,kz).gt.smallnum) then
	464	if (sp_zer.eqv..true.) then ! first non zero value
	465	sp_count_i=sp_count_i+1
	466	sparse_dump_i(sp_count_i)= &
	467	ix+jynumxgrid+kznumxgrid*numygrid
	468	sp_zer=.false.
	469	sp_fact=sp_fact*(-1.)
	470	endif
	471	sp_count_r=sp_count_r+1
[20963b1]	472	if (lparticlecountoutput) then
	473	sparse_dump_r(sp_count_r)= &
	474	sp_fact* &
	475	grid(ix,jy,kz)
	476	else
[8a65cb0]	477	sparse_dump_r(sp_count_r)= &
	478	sp_fact* &
	479	grid(ix,jy,kz)* &
	480	factor3d(ix,jy,kz)/tot_mu(ks,kp)
[20963b1]	481	end if
	482
	483
[8a65cb0]	484	! if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
	485	! + write (,) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
	486	! sparse_dump_u(sp_count_r)=
	487	!+ ,gridsigma(ix,jy,kz,ks,kp,nage)*
	488	!+ factor(ix,jy,kz)/tot_mu(ks,kp)
	489	else ! concentration is zero
	490	sp_zer=.true.
	491	endif
	492	end do
	493	end do
	494	end do
	495	write(unitoutgrid) sp_count_i
	496	write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
	497	write(unitoutgrid) sp_count_r
	498	write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
	499	! write(unitoutgrid) sp_count_u
	500	! write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
	501
	502
	503
	504	endif ! concentration output
	505
	506	! Mixing ratio output
	507	!********************
	508
	509	if ((iout.eq.2).or.(iout.eq.3)) then ! mixing ratio
	510
	511	! Wet deposition
	512	sp_count_i=0
	513	sp_count_r=0
	514	sp_fact=-1.
	515	sp_zer=.true.
	516	if ((ldirect.eq.1).and.(WETDEP)) then
	517	do jy=0,numygrid-1
	518	do ix=0,numxgrid-1
	519	if (wetgrid(ix,jy).gt.smallnum) then
	520	if (sp_zer.eqv..true.) then ! first non zero value
	521	sp_count_i=sp_count_i+1
	522	sparse_dump_i(sp_count_i)= &
	523	ix+jy*numxgrid
	524	sp_zer=.false.
	525	sp_fact=sp_fact*(-1.)
	526	endif
	527	sp_count_r=sp_count_r+1
	528	sparse_dump_r(sp_count_r)= &
	529	sp_fact* &
	530	1.e12*wetgrid(ix,jy)/area(ix,jy)
	531	! sparse_dump_u(sp_count_r)=
	532	! + ,1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
	533	else ! concentration is zero
	534	sp_zer=.true.
	535	endif
	536	end do
	537	end do
	538	else
	539	sp_count_i=0
	540	sp_count_r=0
	541	endif
	542	write(unitoutgridppt) sp_count_i
	543	write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
	544	write(unitoutgridppt) sp_count_r
	545	write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
	546	! write(unitoutgridppt) sp_count_u
	547	! write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
	548
	549
	550	! Dry deposition
	551	sp_count_i=0
	552	sp_count_r=0
	553	sp_fact=-1.
	554	sp_zer=.true.
	555	if ((ldirect.eq.1).and.(DRYDEP)) then
	556	do jy=0,numygrid-1
	557	do ix=0,numxgrid-1
	558	if (drygrid(ix,jy).gt.smallnum) then
	559	if (sp_zer.eqv..true.) then ! first non zero value
	560	sp_count_i=sp_count_i+1
	561	sparse_dump_i(sp_count_i)= &
	562	ix+jy*numxgrid
	563	sp_zer=.false.
	564	sp_fact=sp_fact*(-1)
	565	endif
	566	sp_count_r=sp_count_r+1
	567	sparse_dump_r(sp_count_r)= &
	568	sp_fact* &
	569	1.e12*drygrid(ix,jy)/area(ix,jy)
	570	! sparse_dump_u(sp_count_r)=
	571	! + ,1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
	572	else ! concentration is zero
	573	sp_zer=.true.
	574	endif
	575	end do
	576	end do
	577	else
	578	sp_count_i=0
	579	sp_count_r=0
	580	endif
	581	write(unitoutgridppt) sp_count_i
	582	write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
	583	write(unitoutgridppt) sp_count_r
	584	write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
	585	! write(unitoutgridppt) sp_count_u
	586	! write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
	587
	588
	589	! Mixing ratios
	590	sp_count_i=0
	591	sp_count_r=0
	592	sp_fact=-1.
	593	sp_zer=.true.
	594	do kz=1,numzgrid
	595	do jy=0,numygrid-1
	596	do ix=0,numxgrid-1
	597	if (grid(ix,jy,kz).gt.smallnum) then
	598	if (sp_zer.eqv..true.) then ! first non zero value
	599	sp_count_i=sp_count_i+1
	600	sparse_dump_i(sp_count_i)= &
	601	ix+jynumxgrid+kznumxgrid*numygrid
	602	sp_zer=.false.
	603	sp_fact=sp_fact*(-1.)
	604	endif
	605	sp_count_r=sp_count_r+1
	606	sparse_dump_r(sp_count_r)= &
	607	sp_fact* &
	608	1.e12*grid(ix,jy,kz) &
	609	/volume(ix,jy,kz)/outnum* &
	610	weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
	611	! sparse_dump_u(sp_count_r)=
	612	!+ ,1.e12*gridsigma(ix,jy,kz,ks,kp,nage)/volume(ix,jy,kz)/
	613	!+ outnum*weightair/weightmolar(ks)/
	614	!+ densityoutgrid(ix,jy,kz)
	615	else ! concentration is zero
	616	sp_zer=.true.
	617	endif
	618	end do
	619	end do
	620	end do
	621	write(unitoutgridppt) sp_count_i
	622	write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
	623	write(unitoutgridppt) sp_count_r
	624	write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
	625	! write(unitoutgridppt) sp_count_u
	626	! write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
	627
	628	endif ! output for ppt
	629
	630	end do
	631	end do
	632
	633	close(unitoutgridppt)
	634	close(unitoutgrid)
	635
	636	end do
	637
[f3054ea]	638	! RLT Aug 2017
	639	! Write out conversion factor for dry air
	640	inquire(file=path(2)(1:length(2))//'factor_drygrid',exist=lexist)
	641	if (lexist) then
	642	! open and append
	643	open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid',form='unformatted',&
	644	status='old',action='write',access='append')
	645	else
	646	! create new
	647	open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid',form='unformatted',&
	648	status='new',action='write')
	649	endif
	650	sp_count_i=0
	651	sp_count_r=0
	652	sp_fact=-1.
	653	sp_zer=.true.
	654	do kz=1,numzgrid
	655	do jy=0,numygrid-1
	656	do ix=0,numxgrid-1
	657	if (factor_drygrid(ix,jy,kz).gt.(1.+smallnum).or.factor_drygrid(ix,jy,kz).lt.(1.-smallnum)) then
	658	if (sp_zer.eqv..true.) then ! first value not equal to one
	659	sp_count_i=sp_count_i+1
	660	sparse_dump_i(sp_count_i)= &
	661	ix+jynumxgrid+kznumxgrid*numygrid
	662	sp_zer=.false.
	663	sp_fact=sp_fact*(-1.)
	664	endif
	665	sp_count_r=sp_count_r+1
	666	sparse_dump_r(sp_count_r)= &
	667	sp_fact*factor_drygrid(ix,jy,kz)
	668	else ! factor is one
	669	sp_zer=.true.
	670	endif
	671	end do
	672	end do
	673	end do
	674	write(unitoutfactor) sp_count_i
	675	write(unitoutfactor) (sparse_dump_i(i),i=1,sp_count_i)
	676	write(unitoutfactor) sp_count_r
	677	write(unitoutfactor) (sparse_dump_r(i),i=1,sp_count_r)
	678	close(unitoutfactor)
	679
	680
[8a65cb0]	681	if (gridtotal.gt.0.) gridtotalunc=gridsigmatotal/gridtotal
	682	if (wetgridtotal.gt.0.) wetgridtotalunc=wetgridsigmatotal/ &
	683	wetgridtotal
	684	if (drygridtotal.gt.0.) drygridtotalunc=drygridsigmatotal/ &
	685	drygridtotal
	686
	687	! Dump of receptor concentrations
	688
	689	if (numreceptor.gt.0 .and. (iout.eq.2 .or. iout.eq.3) ) then
	690	write(unitoutreceptppt) itime
	691	do ks=1,nspec
	692	write(unitoutreceptppt) (1.e12creceptor(i,ks)/outnum &
	693	weightair/weightmolar(ks)/densityoutrecept(i),i=1,numreceptor)
	694	end do
	695	endif
	696
	697	! Dump of receptor concentrations
	698
	699	if (numreceptor.gt.0) then
	700	write(unitoutrecept) itime
	701	do ks=1,nspec
	702	write(unitoutrecept) (1.e12*creceptor(i,ks)/outnum, &
	703	i=1,numreceptor)
	704	end do
	705	endif
	706
[f3054ea]	707	! RLT Aug 2017
	708	! Write out conversion factor for dry air
	709	if (numreceptor.gt.0) then
	710	inquire(file=path(2)(1:length(2))//'factor_dryreceptor',exist=lexist)
	711	if (lexist) then
	712	! open and append
	713	open(unitoutfactor,file=path(2)(1:length(2))//'factor_dryreceptor',form='unformatted',&
	714	status='old',action='write',access='append')
	715	else
	716	! create new
	717	open(unitoutfactor,file=path(2)(1:length(2))//'factor_dryreceptor',form='unformatted',&
	718	status='new',action='write')
	719	endif
	720	write(unitoutfactor) itime
	721	write(unitoutfactor) (factor_dryrecept(i),i=1,numreceptor)
	722	close(unitoutfactor)
	723	endif
[8a65cb0]	724
	725	! Reinitialization of grid
	726	!*************************
	727
[20963b1]	728	creceptor(:,:)=0.
	729	gridunc(:,:,:,:,:,:,:)=0.
[8a65cb0]	730
	731	if (mp_measure_time) call mpif_mtime('rootonly',1)
	732
	733	end subroutine concoutput

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

Download in other formats:

Original Format