source: flexpart.git/src/concoutput_mpi.f90 @ fddc6ec

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bugunivie
Last change on this file since fddc6ec was fdc0f03, checked in by Espen Sollum ATMOS <eso@…>, 8 years ago

Code for cloud water should be correct if the total cw is stored in field clwc (old scheme) or in field qc (new scheme). Minor edits in some files.

  • Property mode set to 100644
File size: 23.9 KB
Line 
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
22subroutine 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
70
71  implicit none
72
73  real(kind=dp) :: jul
74  integer :: itime,i,ix,jy,kz,ks,kp,l,iix,jjy,kzz,nage,jjjjmmdd,ihmmss
75  integer :: sp_count_i,sp_count_r
76  real :: sp_fact
77  real :: outnum,densityoutrecept(maxreceptor),xl,yl
78
79!real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
80!    +grid(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,maxpointspec_act,
81!    +    maxageclass)
82!real wetgrid(0:numxgrid-1,0:numygrid-1,maxspec,maxpointspec_act,
83!    +       maxageclass)
84!real drygrid(0:numxgrid-1,0:numygrid-1,maxspec,
85!    +       maxpointspec_act,maxageclass)
86!real gridsigma(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,
87!    +       maxpointspec_act,maxageclass),
88!    +     drygridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
89!    +     maxpointspec_act,maxageclass),
90!    +     wetgridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
91!    +     maxpointspec_act,maxageclass)
92!real factor(0:numxgrid-1,0:numygrid-1,numzgrid)
93!real sparse_dump_r(numxgrid*numygrid*numzgrid)
94!integer sparse_dump_i(numxgrid*numygrid*numzgrid)
95
96!real sparse_dump_u(numxgrid*numygrid*numzgrid)
97  real :: auxgrid(nclassunc),gridtotal,gridsigmatotal,gridtotalunc
98  real :: wetgridtotal,wetgridsigmatotal,wetgridtotalunc
99  real :: drygridtotal,drygridsigmatotal,drygridtotalunc
100  real :: halfheight,dz,dz1,dz2,tot_mu(maxspec,maxpointspec_act)
101  real,parameter :: smallnum = tiny(0.0) ! smallest number that can be handled
102  real,parameter :: weightair=28.97
103  logical :: sp_zer
104  LOGICAL,save :: init=.true.
105  character :: adate*8,atime*6
106  character(len=3) :: anspec
107  integer :: mind
108! mind        eso: added to ensure identical results between 2&3-fields versions
109  character(LEN=8),save :: file_stat='REPLACE'
110  logical :: ldates_file
111  integer :: ierr
112  character(LEN=100) :: dates_char
113!  character :: dates_char
114
115! Measure execution time
116  if (mp_measure_time) call mpif_mtime('rootonly',0)
117
118! Determine current calendar date, needed for the file name
119!**********************************************************
120
121  jul=bdate+real(itime,kind=dp)/86400._dp
122  call caldate(jul,jjjjmmdd,ihmmss)
123  write(adate,'(i8.8)') jjjjmmdd
124  write(atime,'(i6.6)') ihmmss
125
126! Overwrite existing dates file on first call, later append to it
127! This fixes a bug where the dates file kept growing across multiple runs
128
129! If 'dates' file exists, make a backup
130  inquire(file=path(2)(1:length(2))//'dates', exist=ldates_file)
131  if (ldates_file.and.init) then
132    open(unit=unitdates, file=path(2)(1:length(2))//'dates',form='formatted', &
133         &access='sequential', status='old', action='read', iostat=ierr)
134    open(unit=unittmp, file=path(2)(1:length(2))//'dates.bak', access='sequential', &
135         &status='replace', action='write', form='formatted', iostat=ierr)
136    do while (.true.)
137      read(unitdates, '(a)', iostat=ierr) dates_char
138      if (ierr.ne.0) exit
139      write(unit=unittmp, fmt='(a)', iostat=ierr, advance='yes') trim(dates_char)
140!      if (ierr.ne.0) write(*,*) "Write error, ", ierr
141    end do
142    close(unit=unitdates)
143    close(unit=unittmp)
144  end if
145
146  open(unitdates,file=path(2)(1:length(2))//'dates', ACCESS='APPEND', STATUS=file_stat)
147  write(unitdates,'(a)') adate//atime
148  close(unitdates) 
149
150  ! Overwrite existing dates file on first call, later append to it
151  ! This fixes a bug where the dates file kept growing across multiple runs
152  IF (init) THEN
153    file_stat='OLD'
154    init=.false.
155  END IF
156
157
158! For forward simulations, output fields have dimension MAXSPEC,
159! for backward simulations, output fields have dimension MAXPOINT.
160! Thus, make loops either about nspec, or about numpoint
161!*****************************************************************
162
163
164  if (ldirect.eq.1) then
165    do ks=1,nspec
166      do kp=1,maxpointspec_act
167        tot_mu(ks,kp)=1
168      end do
169    end do
170  else
171    do ks=1,nspec
172      do kp=1,maxpointspec_act
173        tot_mu(ks,kp)=xmass(kp,ks)
174      end do
175    end do
176  endif
177
178
179!*******************************************************************
180! Compute air density: sufficiently accurate to take it
181! from coarse grid at some time
182! Determine center altitude of output layer, and interpolate density
183! data to that altitude
184!*******************************************************************
185
186  mind=memind(2)
187  do kz=1,numzgrid
188    if (kz.eq.1) then
189      halfheight=outheight(1)/2.
190    else
191      halfheight=(outheight(kz)+outheight(kz-1))/2.
192    endif
193    do kzz=2,nz
194      if ((height(kzz-1).lt.halfheight).and. &
195           (height(kzz).gt.halfheight)) goto 46
196    end do
19746  kzz=max(min(kzz,nz),2)
198    dz1=halfheight-height(kzz-1)
199    dz2=height(kzz)-halfheight
200    dz=dz1+dz2
201    do jy=0,numygrid-1
202      do ix=0,numxgrid-1
203        xl=outlon0+real(ix)*dxout
204        yl=outlat0+real(jy)*dyout
205        xl=(xl-xlon0)/dx
206        yl=(yl-ylat0)/dy !v9.1.1
207        iix=max(min(nint(xl),nxmin1),0)
208        jjy=max(min(nint(yl),nymin1),0)
209        ! densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
210        !      rho(iix,jjy,kzz-1,2)*dz2)/dz
211        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,mind)*dz1+ &
212             rho(iix,jjy,kzz-1,mind)*dz2)/dz
213      end do
214    end do
215  end do
216
217  do i=1,numreceptor
218    xl=xreceptor(i)
219    yl=yreceptor(i)
220    iix=max(min(nint(xl),nxmin1),0)
221    jjy=max(min(nint(yl),nymin1),0)
222    !densityoutrecept(i)=rho(iix,jjy,1,2)
223    densityoutrecept(i)=rho(iix,jjy,1,mind)
224  end do
225
226
227! Output is different for forward and backward simulations
228  do kz=1,numzgrid
229    do jy=0,numygrid-1
230      do ix=0,numxgrid-1
231        if (ldirect.eq.1) then
232          factor3d(ix,jy,kz)=1.e12/volume(ix,jy,kz)/outnum
233        else
234          factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
235        endif
236      end do
237    end do
238  end do
239
240!*********************************************************************
241! Determine the standard deviation of the mean concentration or mixing
242! ratio (uncertainty of the output) and the dry and wet deposition
243!*********************************************************************
244
245  gridtotal=0.
246  gridsigmatotal=0.
247  gridtotalunc=0.
248  wetgridtotal=0.
249  wetgridsigmatotal=0.
250  wetgridtotalunc=0.
251  drygridtotal=0.
252  drygridsigmatotal=0.
253  drygridtotalunc=0.
254
255  do ks=1,nspec
256
257    write(anspec,'(i3.3)') ks
258    if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
259      if (ldirect.eq.1) then
260        open(unitoutgrid,file=path(2)(1:length(2))//'grid_conc_'//adate// &
261             atime//'_'//anspec,form='unformatted')
262      else
263        open(unitoutgrid,file=path(2)(1:length(2))//'grid_time_'//adate// &
264             atime//'_'//anspec,form='unformatted')
265      endif
266      write(unitoutgrid) itime
267    endif
268
269    if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
270      open(unitoutgridppt,file=path(2)(1:length(2))//'grid_pptv_'//adate// &
271           atime//'_'//anspec,form='unformatted')
272
273      write(unitoutgridppt) itime
274    endif
275
276    do kp=1,maxpointspec_act
277      do nage=1,nageclass
278
279        do jy=0,numygrid-1
280          do ix=0,numxgrid-1
281
282! WET DEPOSITION
283            if ((WETDEP).and.(ldirect.gt.0)) then
284              do l=1,nclassunc
285                auxgrid(l)=wetgridunc0(ix,jy,ks,kp,l,nage)
286              end do
287              call mean(auxgrid,wetgrid(ix,jy), &
288                   wetgridsigma(ix,jy),nclassunc)
289! Multiply by number of classes to get total concentration
290              wetgrid(ix,jy)=wetgrid(ix,jy) &
291                   *nclassunc
292              wetgridtotal=wetgridtotal+wetgrid(ix,jy)
293! Calculate standard deviation of the mean
294              wetgridsigma(ix,jy)= &
295                   wetgridsigma(ix,jy)* &
296                   sqrt(real(nclassunc))
297              wetgridsigmatotal=wetgridsigmatotal+ &
298                   wetgridsigma(ix,jy)
299            endif
300
301! DRY DEPOSITION
302            if ((DRYDEP).and.(ldirect.gt.0)) then
303              do l=1,nclassunc
304                auxgrid(l)=drygridunc0(ix,jy,ks,kp,l,nage)
305              end do
306              call mean(auxgrid,drygrid(ix,jy), &
307                   drygridsigma(ix,jy),nclassunc)
308! Multiply by number of classes to get total concentration
309              drygrid(ix,jy)=drygrid(ix,jy)* &
310                   nclassunc
311              drygridtotal=drygridtotal+drygrid(ix,jy)
312! Calculate standard deviation of the mean
313              drygridsigma(ix,jy)= &
314                   drygridsigma(ix,jy)* &
315                   sqrt(real(nclassunc))
316              drygridsigmatotal=drygridsigmatotal+ &
317                   drygridsigma(ix,jy)
318            endif
319
320! CONCENTRATION OR MIXING RATIO
321            do kz=1,numzgrid
322              do l=1,nclassunc
323                auxgrid(l)=gridunc(ix,jy,kz,ks,kp,l,nage)
324              end do
325              call mean(auxgrid,grid(ix,jy,kz), &
326                   gridsigma(ix,jy,kz),nclassunc)
327! Multiply by number of classes to get total concentration
328              grid(ix,jy,kz)= &
329                   grid(ix,jy,kz)*nclassunc
330              gridtotal=gridtotal+grid(ix,jy,kz)
331! Calculate standard deviation of the mean
332              gridsigma(ix,jy,kz)= &
333                   gridsigma(ix,jy,kz)* &
334                   sqrt(real(nclassunc))
335              gridsigmatotal=gridsigmatotal+ &
336                   gridsigma(ix,jy,kz)
337            end do
338          end do
339        end do
340
341
342
343
344!*******************************************************************
345! Generate output: may be in concentration (ng/m3) or in mixing
346! ratio (ppt) or both
347! Output the position and the values alternated multiplied by
348! 1 or -1, first line is number of values, number of positions
349! For backward simulations, the unit is seconds, stored in grid_time
350!*******************************************************************
351
352! Concentration output
353!*********************
354        if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
355
356! Wet deposition
357          sp_count_i=0
358          sp_count_r=0
359          sp_fact=-1.
360          sp_zer=.true.
361          if ((ldirect.eq.1).and.(WETDEP)) then
362            do jy=0,numygrid-1
363              do ix=0,numxgrid-1
364!oncentraion greater zero
365                if (wetgrid(ix,jy).gt.smallnum) then
366                  if (sp_zer.eqv..true.) then ! first non zero value
367                    sp_count_i=sp_count_i+1
368                    sparse_dump_i(sp_count_i)=ix+jy*numxgrid
369                    sp_zer=.false.
370                    sp_fact=sp_fact*(-1.)
371                  endif
372                  sp_count_r=sp_count_r+1
373                  sparse_dump_r(sp_count_r)= &
374                       sp_fact*1.e12*wetgrid(ix,jy)/area(ix,jy)
375!                sparse_dump_u(sp_count_r)=
376!+                1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
377                else ! concentration is zero
378                  sp_zer=.true.
379                endif
380              end do
381            end do
382          else
383            sp_count_i=0
384            sp_count_r=0
385          endif
386          write(unitoutgrid) sp_count_i
387          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
388          write(unitoutgrid) sp_count_r
389          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
390!       write(unitoutgrid) sp_count_u
391!       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
392
393! Dry deposition
394          sp_count_i=0
395          sp_count_r=0
396          sp_fact=-1.
397          sp_zer=.true.
398          if ((ldirect.eq.1).and.(DRYDEP)) then
399            do jy=0,numygrid-1
400              do ix=0,numxgrid-1
401                if (drygrid(ix,jy).gt.smallnum) then
402                  if (sp_zer.eqv..true.) then ! first non zero value
403                    sp_count_i=sp_count_i+1
404                    sparse_dump_i(sp_count_i)=ix+jy*numxgrid
405                    sp_zer=.false.
406                    sp_fact=sp_fact*(-1.)
407                  endif
408                  sp_count_r=sp_count_r+1
409                  sparse_dump_r(sp_count_r)= &
410                       sp_fact* &
411                       1.e12*drygrid(ix,jy)/area(ix,jy)
412!                sparse_dump_u(sp_count_r)=
413!+                1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
414                else ! concentration is zero
415                  sp_zer=.true.
416                endif
417              end do
418            end do
419          else
420            sp_count_i=0
421            sp_count_r=0
422          endif
423          write(unitoutgrid) sp_count_i
424          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
425          write(unitoutgrid) sp_count_r
426          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
427!       write(*,*) sp_count_u
428!       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
429
430
431
432! Concentrations
433          sp_count_i=0
434          sp_count_r=0
435          sp_fact=-1.
436          sp_zer=.true.
437          do kz=1,numzgrid
438            do jy=0,numygrid-1
439              do ix=0,numxgrid-1
440                if (grid(ix,jy,kz).gt.smallnum) then
441                  if (sp_zer.eqv..true.) then ! first non zero value
442                    sp_count_i=sp_count_i+1
443                    sparse_dump_i(sp_count_i)= &
444                         ix+jy*numxgrid+kz*numxgrid*numygrid
445                    sp_zer=.false.
446                    sp_fact=sp_fact*(-1.)
447                  endif
448                  sp_count_r=sp_count_r+1
449                  sparse_dump_r(sp_count_r)= &
450                       sp_fact* &
451                       grid(ix,jy,kz)* &
452                       factor3d(ix,jy,kz)/tot_mu(ks,kp)
453!                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
454!    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
455!                sparse_dump_u(sp_count_r)=
456!+               ,gridsigma(ix,jy,kz,ks,kp,nage)*
457!+               factor(ix,jy,kz)/tot_mu(ks,kp)
458                else ! concentration is zero
459                  sp_zer=.true.
460                endif
461              end do
462            end do
463          end do
464          write(unitoutgrid) sp_count_i
465          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
466          write(unitoutgrid) sp_count_r
467          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
468!       write(unitoutgrid) sp_count_u
469!       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
470
471
472
473        endif !  concentration output
474
475! Mixing ratio output
476!********************
477
478        if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
479
480! Wet deposition
481          sp_count_i=0
482          sp_count_r=0
483          sp_fact=-1.
484          sp_zer=.true.
485          if ((ldirect.eq.1).and.(WETDEP)) then
486            do jy=0,numygrid-1
487              do ix=0,numxgrid-1
488                if (wetgrid(ix,jy).gt.smallnum) then
489                  if (sp_zer.eqv..true.) then ! first non zero value
490                    sp_count_i=sp_count_i+1
491                    sparse_dump_i(sp_count_i)= &
492                         ix+jy*numxgrid
493                    sp_zer=.false.
494                    sp_fact=sp_fact*(-1.)
495                  endif
496                  sp_count_r=sp_count_r+1
497                  sparse_dump_r(sp_count_r)= &
498                       sp_fact* &
499                       1.e12*wetgrid(ix,jy)/area(ix,jy)
500!                sparse_dump_u(sp_count_r)=
501!    +            ,1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
502                else ! concentration is zero
503                  sp_zer=.true.
504                endif
505              end do
506            end do
507          else
508            sp_count_i=0
509            sp_count_r=0
510          endif
511          write(unitoutgridppt) sp_count_i
512          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
513          write(unitoutgridppt) sp_count_r
514          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
515!       write(unitoutgridppt) sp_count_u
516!       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
517
518
519! Dry deposition
520          sp_count_i=0
521          sp_count_r=0
522          sp_fact=-1.
523          sp_zer=.true.
524          if ((ldirect.eq.1).and.(DRYDEP)) then
525            do jy=0,numygrid-1
526              do ix=0,numxgrid-1
527                if (drygrid(ix,jy).gt.smallnum) then
528                  if (sp_zer.eqv..true.) then ! first non zero value
529                    sp_count_i=sp_count_i+1
530                    sparse_dump_i(sp_count_i)= &
531                         ix+jy*numxgrid
532                    sp_zer=.false.
533                    sp_fact=sp_fact*(-1)
534                  endif
535                  sp_count_r=sp_count_r+1
536                  sparse_dump_r(sp_count_r)= &
537                       sp_fact* &
538                       1.e12*drygrid(ix,jy)/area(ix,jy)
539!                sparse_dump_u(sp_count_r)=
540!    +            ,1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
541                else ! concentration is zero
542                  sp_zer=.true.
543                endif
544              end do
545            end do
546          else
547            sp_count_i=0
548            sp_count_r=0
549          endif
550          write(unitoutgridppt) sp_count_i
551          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
552          write(unitoutgridppt) sp_count_r
553          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
554!       write(unitoutgridppt) sp_count_u
555!       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
556
557
558! Mixing ratios
559          sp_count_i=0
560          sp_count_r=0
561          sp_fact=-1.
562          sp_zer=.true.
563          do kz=1,numzgrid
564            do jy=0,numygrid-1
565              do ix=0,numxgrid-1
566                if (grid(ix,jy,kz).gt.smallnum) then
567                  if (sp_zer.eqv..true.) then ! first non zero value
568                    sp_count_i=sp_count_i+1
569                    sparse_dump_i(sp_count_i)= &
570                         ix+jy*numxgrid+kz*numxgrid*numygrid
571                    sp_zer=.false.
572                    sp_fact=sp_fact*(-1.)
573                  endif
574                  sp_count_r=sp_count_r+1
575                  sparse_dump_r(sp_count_r)= &
576                       sp_fact* &
577                       1.e12*grid(ix,jy,kz) &
578                       /volume(ix,jy,kz)/outnum* &
579                       weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
580!                sparse_dump_u(sp_count_r)=
581!+              ,1.e12*gridsigma(ix,jy,kz,ks,kp,nage)/volume(ix,jy,kz)/
582!+              outnum*weightair/weightmolar(ks)/
583!+              densityoutgrid(ix,jy,kz)
584                else ! concentration is zero
585                  sp_zer=.true.
586                endif
587              end do
588            end do
589          end do
590          write(unitoutgridppt) sp_count_i
591          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
592          write(unitoutgridppt) sp_count_r
593          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
594!       write(unitoutgridppt) sp_count_u
595!       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
596
597        endif ! output for ppt
598
599      end do
600    end do
601
602    close(unitoutgridppt)
603    close(unitoutgrid)
604
605  end do
606
607  if (gridtotal.gt.0.) gridtotalunc=gridsigmatotal/gridtotal
608  if (wetgridtotal.gt.0.) wetgridtotalunc=wetgridsigmatotal/ &
609       wetgridtotal
610  if (drygridtotal.gt.0.) drygridtotalunc=drygridsigmatotal/ &
611       drygridtotal
612
613! Dump of receptor concentrations
614
615  if (numreceptor.gt.0 .and. (iout.eq.2 .or. iout.eq.3)  ) then
616    write(unitoutreceptppt) itime
617    do ks=1,nspec
618      write(unitoutreceptppt) (1.e12*creceptor(i,ks)/outnum* &
619           weightair/weightmolar(ks)/densityoutrecept(i),i=1,numreceptor)
620    end do
621  endif
622
623! Dump of receptor concentrations
624
625  if (numreceptor.gt.0) then
626    write(unitoutrecept) itime
627    do ks=1,nspec
628      write(unitoutrecept) (1.e12*creceptor(i,ks)/outnum, &
629           i=1,numreceptor)
630    end do
631  endif
632
633
634
635! Reinitialization of grid
636!*************************
637
638  do ks=1,nspec
639    do kp=1,maxpointspec_act
640      do i=1,numreceptor
641        creceptor(i,ks)=0.
642      end do
643      do jy=0,numygrid-1
644        do ix=0,numxgrid-1
645          do l=1,nclassunc
646            do nage=1,nageclass
647              do kz=1,numzgrid
648                gridunc(ix,jy,kz,ks,kp,l,nage)=0.
649              end do
650            end do
651          end do
652        end do
653      end do
654    end do
655  end do
656
657  if (mp_measure_time) call mpif_mtime('rootonly',1)
658 
659end subroutine concoutput
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG