source: flexpart.git/src/concoutput_surf_mpi.f90 @ 38b7917

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

Added option to use double precision for calculating the deposition fields

  • Property mode set to 100644
File size: 23.6 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_surf(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  !                                                                            *
57  ! Variables:                                                                 *
58  ! outnum          number of samples                                          *
59  ! ncells          number of cells with non-zero concentrations               *
60  ! sparse          .true. if in sparse matrix format, else .false.            *
61  ! tot_mu          1 for forward, initial mass mixing ration for backw. runs  *
62  !                                                                            *
63  !*****************************************************************************
64
65  use unc_mod
66  use point_mod
67  use outg_mod
68  use par_mod
69  use com_mod
70  use mpi_mod
71  use mean_mod
72
73  implicit none
74
75  real(kind=dp) :: jul
76  integer :: itime,i,ix,jy,kz,ks,kp,l,iix,jjy,kzz,nage,jjjjmmdd,ihmmss
77  integer :: sp_count_i,sp_count_r
78  real :: sp_fact
79  real :: outnum,densityoutrecept(maxreceptor),xl,yl
80
81  !real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
82  !    +grid(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,maxpointspec_act,
83  !    +    maxageclass)
84  !real wetgrid(0:numxgrid-1,0:numygrid-1,maxspec,maxpointspec_act,
85  !    +       maxageclass)
86  !real drygrid(0:numxgrid-1,0:numygrid-1,maxspec,
87  !    +       maxpointspec_act,maxageclass)
88  !real gridsigma(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,
89  !    +       maxpointspec_act,maxageclass),
90  !    +     drygridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
91  !    +     maxpointspec_act,maxageclass),
92  !    +     wetgridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
93  !    +     maxpointspec_act,maxageclass)
94  !real factor(0:numxgrid-1,0:numygrid-1,numzgrid)
95  !real sparse_dump_r(numxgrid*numygrid*numzgrid)
96  !integer sparse_dump_i(numxgrid*numygrid*numzgrid)
97
98  !real sparse_dump_u(numxgrid*numygrid*numzgrid)
99  real(dep_prec) :: auxgrid(nclassunc)
100  real(sp) :: gridtotal,gridsigmatotal,gridtotalunc
101  real(dep_prec) :: wetgridtotal,wetgridsigmatotal,wetgridtotalunc
102  real(dep_prec) :: drygridtotal,drygridsigmatotal,drygridtotalunc
103  real :: halfheight,dz,dz1,dz2,tot_mu(maxspec,maxpointspec_act)
104  real,parameter :: smallnum = tiny(0.0) ! smallest number that can be handled
105  real,parameter :: weightair=28.97
106  logical :: sp_zer
107  character :: adate*8,atime*6
108  character(len=3) :: anspec
109  integer :: mind
110! mind        eso:added to get consistent results between 2&3-fields versions
111
112! Measure execution time
113  if (mp_measure_time) call mpif_mtime('rootonly',0)
114
115
116  ! Determine current calendar date, needed for the file name
117  !**********************************************************
118
119  jul=bdate+real(itime,kind=dp)/86400._dp
120  call caldate(jul,jjjjmmdd,ihmmss)
121  write(adate,'(i8.8)') jjjjmmdd
122  write(atime,'(i6.6)') ihmmss
123  !write(unitdates,'(a)') adate//atime
124
125    open(unitdates,file=path(2)(1:length(2))//'dates', ACCESS='APPEND')
126      write(unitdates,'(a)') adate//atime
127    close(unitdates)
128
129  ! For forward simulations, output fields have dimension MAXSPEC,
130  ! for backward simulations, output fields have dimension MAXPOINT.
131  ! Thus, make loops either about nspec, or about numpoint
132  !*****************************************************************
133
134
135  if (ldirect.eq.1) then
136    do ks=1,nspec
137      do kp=1,maxpointspec_act
138        tot_mu(ks,kp)=1
139      end do
140    end do
141  else
142    do ks=1,nspec
143      do kp=1,maxpointspec_act
144        tot_mu(ks,kp)=xmass(kp,ks)
145      end do
146    end do
147  endif
148
149
150  if (verbosity.eq.1) then
151     print*,'concoutput_surf 2'
152     CALL SYSTEM_CLOCK(count_clock)
153     WRITE(*,*) 'SYSTEM_CLOCK',count_clock - count_clock0   
154  endif
155
156  !*******************************************************************
157  ! Compute air density: sufficiently accurate to take it
158  ! from coarse grid at some time
159  ! Determine center altitude of output layer, and interpolate density
160  ! data to that altitude
161  !*******************************************************************
162
163  mind=memind(2)
164  do kz=1,numzgrid
165    if (kz.eq.1) then
166      halfheight=outheight(1)/2.
167    else
168      halfheight=(outheight(kz)+outheight(kz-1))/2.
169    endif
170    do kzz=2,nz
171      if ((height(kzz-1).lt.halfheight).and. &
172           (height(kzz).gt.halfheight)) goto 46
173    end do
17446   kzz=max(min(kzz,nz),2)
175    dz1=halfheight-height(kzz-1)
176    dz2=height(kzz)-halfheight
177    dz=dz1+dz2
178    do jy=0,numygrid-1
179      do ix=0,numxgrid-1
180        xl=outlon0+real(ix)*dxout
181        yl=outlat0+real(jy)*dyout
182        xl=(xl-xlon0)/dx
183        yl=(yl-ylat0)/dy
184        iix=max(min(nint(xl),nxmin1),0)
185        jjy=max(min(nint(yl),nymin1),0)
186        ! densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
187        !      rho(iix,jjy,kzz-1,2)*dz2)/dz
188        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,mind)*dz1+ &
189             rho(iix,jjy,kzz-1,mind)*dz2)/dz
190      end do
191    end do
192  end do
193
194  do i=1,numreceptor
195    xl=xreceptor(i)
196    yl=yreceptor(i)
197    iix=max(min(nint(xl),nxmin1),0)
198    jjy=max(min(nint(yl),nymin1),0)
199    !densityoutrecept(i)=rho(iix,jjy,1,2)
200    densityoutrecept(i)=rho(iix,jjy,1,mind)
201  end do
202
203
204  ! Output is different for forward and backward simulations
205    do kz=1,numzgrid
206      do jy=0,numygrid-1
207        do ix=0,numxgrid-1
208          if (ldirect.eq.1) then
209            factor3d(ix,jy,kz)=1.e12/volume(ix,jy,kz)/outnum
210          else
211            factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
212          endif
213        end do
214      end do
215    end do
216
217  !*********************************************************************
218  ! Determine the standard deviation of the mean concentration or mixing
219  ! ratio (uncertainty of the output) and the dry and wet deposition
220  !*********************************************************************
221
222  if (verbosity.eq.1) then
223     print*,'concoutput_surf 3 (sd)'
224     CALL SYSTEM_CLOCK(count_clock)
225     WRITE(*,*) 'SYSTEM_CLOCK',count_clock - count_clock0   
226  endif
227  gridtotal=0.
228  gridsigmatotal=0.
229  gridtotalunc=0.
230  wetgridtotal=0.
231  wetgridsigmatotal=0.
232  wetgridtotalunc=0.
233  drygridtotal=0.
234  drygridsigmatotal=0.
235  drygridtotalunc=0.
236
237  do ks=1,nspec
238
239  write(anspec,'(i3.3)') ks
240  if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
241    if (ldirect.eq.1) then
242      open(unitoutgrid,file=path(2)(1:length(2))//'grid_conc_'//adate// &
243           atime//'_'//anspec,form='unformatted')
244    else
245      open(unitoutgrid,file=path(2)(1:length(2))//'grid_time_'//adate// &
246           atime//'_'//anspec,form='unformatted')
247    endif
248     write(unitoutgrid) itime
249   endif
250
251  if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
252   open(unitoutgridppt,file=path(2)(1:length(2))//'grid_pptv_'//adate// &
253        atime//'_'//anspec,form='unformatted')
254
255    write(unitoutgridppt) itime
256  endif
257
258  do kp=1,maxpointspec_act
259  do nage=1,nageclass
260
261    do jy=0,numygrid-1
262      do ix=0,numxgrid-1
263
264  ! WET DEPOSITION
265        if ((WETDEP).and.(ldirect.gt.0)) then
266            do l=1,nclassunc
267              auxgrid(l)=wetgridunc0(ix,jy,ks,kp,l,nage)
268            end do
269            call mean(auxgrid,wetgrid(ix,jy), &
270                 wetgridsigma(ix,jy),nclassunc)
271  ! Multiply by number of classes to get total concentration
272            wetgrid(ix,jy)=wetgrid(ix,jy) &
273                 *nclassunc
274            wetgridtotal=wetgridtotal+wetgrid(ix,jy)
275  ! Calculate standard deviation of the mean
276            wetgridsigma(ix,jy)= &
277                 wetgridsigma(ix,jy)* &
278                 sqrt(real(nclassunc))
279            wetgridsigmatotal=wetgridsigmatotal+ &
280                 wetgridsigma(ix,jy)
281        endif
282
283  ! DRY DEPOSITION
284        if ((DRYDEP).and.(ldirect.gt.0)) then
285            do l=1,nclassunc
286              auxgrid(l)=drygridunc0(ix,jy,ks,kp,l,nage)
287            end do
288            call mean(auxgrid,drygrid(ix,jy), &
289                 drygridsigma(ix,jy),nclassunc)
290  ! Multiply by number of classes to get total concentration
291            drygrid(ix,jy)=drygrid(ix,jy)* &
292                 nclassunc
293            drygridtotal=drygridtotal+drygrid(ix,jy)
294  ! Calculate standard deviation of the mean
295            drygridsigma(ix,jy)= &
296                 drygridsigma(ix,jy)* &
297                 sqrt(real(nclassunc))
298125         drygridsigmatotal=drygridsigmatotal+ &
299                 drygridsigma(ix,jy)
300        endif
301
302  ! CONCENTRATION OR MIXING RATIO
303        do kz=1,numzgrid
304            do l=1,nclassunc
305              auxgrid(l)=gridunc(ix,jy,kz,ks,kp,l,nage)
306            end do
307            call mean(auxgrid,grid(ix,jy,kz), &
308                 gridsigma(ix,jy,kz),nclassunc)
309  ! Multiply by number of classes to get total concentration
310            grid(ix,jy,kz)= &
311                 grid(ix,jy,kz)*nclassunc
312            gridtotal=gridtotal+grid(ix,jy,kz)
313  ! Calculate standard deviation of the mean
314            gridsigma(ix,jy,kz)= &
315                 gridsigma(ix,jy,kz)* &
316                 sqrt(real(nclassunc))
317            gridsigmatotal=gridsigmatotal+ &
318                 gridsigma(ix,jy,kz)
319        end do
320      end do
321    end do
322
323
324  !*******************************************************************
325  ! Generate output: may be in concentration (ng/m3) or in mixing
326  ! ratio (ppt) or both
327  ! Output the position and the values alternated multiplied by
328  ! 1 or -1, first line is number of values, number of positions
329  ! For backward simulations, the unit is seconds, stored in grid_time
330  !*******************************************************************
331
332  if (verbosity.eq.1) then
333     print*,'concoutput_surf 4 (output)'
334     CALL SYSTEM_CLOCK(count_clock)
335     WRITE(*,*) 'SYSTEM_CLOCK',count_clock - count_clock0   
336  endif
337
338  ! Concentration output
339  !*********************
340
341  if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
342
343  if (verbosity.eq.1) then
344     print*,'concoutput_surf (Wet deposition)'
345     CALL SYSTEM_CLOCK(count_clock)
346     WRITE(*,*) 'SYSTEM_CLOCK',count_clock - count_clock0   
347  endif
348
349  ! Wet deposition
350         sp_count_i=0
351         sp_count_r=0
352         sp_fact=-1.
353         sp_zer=.true.
354         if ((ldirect.eq.1).and.(WETDEP)) then
355         do jy=0,numygrid-1
356            do ix=0,numxgrid-1
357  ! concentraion greater zero
358              if (wetgrid(ix,jy).gt.smallnum) then
359                 if (sp_zer.eqv..true.) then ! first non zero value
360                    sp_count_i=sp_count_i+1
361                    sparse_dump_i(sp_count_i)=ix+jy*numxgrid
362                    sp_zer=.false.
363                    sp_fact=sp_fact*(-1.)
364                 endif
365                 sp_count_r=sp_count_r+1
366                 sparse_dump_r(sp_count_r)= &
367                      sp_fact*1.e12*wetgrid(ix,jy)/area(ix,jy)
368                 sparse_dump_u(sp_count_r)= &
369                      1.e12*wetgridsigma(ix,jy)/area(ix,jy)
370              else ! concentration is zero
371                  sp_zer=.true.
372              endif
373            end do
374         end do
375         else
376            sp_count_i=0
377            sp_count_r=0
378         endif
379         write(unitoutgrid) sp_count_i
380         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
381         write(unitoutgrid) sp_count_r
382         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
383!         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
384
385  if (verbosity.eq.1) then
386     print*,'concoutput_surf (Dry deposition)'
387     CALL SYSTEM_CLOCK(count_clock)
388     WRITE(*,*) 'SYSTEM_CLOCK',count_clock - count_clock0   
389  endif
390  ! Dry deposition
391         sp_count_i=0
392         sp_count_r=0
393         sp_fact=-1.
394         sp_zer=.true.
395         if ((ldirect.eq.1).and.(DRYDEP)) then
396          do jy=0,numygrid-1
397            do ix=0,numxgrid-1
398              if (drygrid(ix,jy).gt.smallnum) then
399                 if (sp_zer.eqv..true.) then ! first non zero value
400                    sp_count_i=sp_count_i+1
401                    sparse_dump_i(sp_count_i)=ix+jy*numxgrid
402                    sp_zer=.false.
403                    sp_fact=sp_fact*(-1.)
404                 endif
405                 sp_count_r=sp_count_r+1
406                 sparse_dump_r(sp_count_r)= &
407                      sp_fact* &
408                      1.e12*drygrid(ix,jy)/area(ix,jy)
409                  sparse_dump_u(sp_count_r)= &
410                      1.e12*drygridsigma(ix,jy)/area(ix,jy)
411              else ! concentration is zero
412                  sp_zer=.true.
413              endif
414            end do
415          end do
416         else
417            sp_count_i=0
418            sp_count_r=0
419         endif
420         write(unitoutgrid) sp_count_i
421         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
422         write(unitoutgrid) sp_count_r
423         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
424!         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
425
426  if (verbosity.eq.1) then
427     print*,'concoutput_surf (Concentrations)'
428     CALL SYSTEM_CLOCK(count_clock)
429     WRITE(*,*) 'SYSTEM_CLOCK',count_clock - count_clock0   
430  endif
431
432  ! Concentrations
433
434  ! surf_only write only 1st layer
435
436         sp_count_i=0
437         sp_count_r=0
438         sp_fact=-1.
439         sp_zer=.true.
440          do kz=1,1
441            do jy=0,numygrid-1
442              do ix=0,numxgrid-1
443                if (grid(ix,jy,kz).gt.smallnum) then
444                  if (sp_zer.eqv..true.) then ! first non zero value
445                    sp_count_i=sp_count_i+1
446                    sparse_dump_i(sp_count_i)= &
447                         ix+jy*numxgrid+kz*numxgrid*numygrid
448                    sp_zer=.false.
449                    sp_fact=sp_fact*(-1.)
450                  endif
451                   sp_count_r=sp_count_r+1
452                   sparse_dump_r(sp_count_r)= &
453                        sp_fact* &
454                        grid(ix,jy,kz)* &
455                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
456  !                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
457  !    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
458                   sparse_dump_u(sp_count_r)= &
459                        gridsigma(ix,jy,kz)* &
460                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
461                else ! concentration is zero
462                  sp_zer=.true.
463                endif
464             end do
465           end do
466         end do
467         write(unitoutgrid) sp_count_i
468         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
469         write(unitoutgrid) sp_count_r
470         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
471!         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
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)/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) (sparse_dump_u(i),i=1,sp_count_r)
516
517
518  ! Dry deposition
519         sp_count_i=0
520         sp_count_r=0
521         sp_fact=-1.
522         sp_zer=.true.
523         if ((ldirect.eq.1).and.(DRYDEP)) then
524          do jy=0,numygrid-1
525            do ix=0,numxgrid-1
526                if (drygrid(ix,jy).gt.smallnum) then
527                  if (sp_zer.eqv..true.) then ! first non zero value
528                    sp_count_i=sp_count_i+1
529                    sparse_dump_i(sp_count_i)= &
530                         ix+jy*numxgrid
531                    sp_zer=.false.
532                    sp_fact=sp_fact*(-1)
533                 endif
534                 sp_count_r=sp_count_r+1
535                 sparse_dump_r(sp_count_r)= &
536                      sp_fact* &
537                      1.e12*drygrid(ix,jy)/area(ix,jy)
538                 sparse_dump_u(sp_count_r)= &
539                      1.e12*drygridsigma(ix,jy)/area(ix,jy)
540              else ! concentration is zero
541                  sp_zer=.true.
542              endif
543            end do
544          end do
545         else
546           sp_count_i=0
547           sp_count_r=0
548         endif
549         write(unitoutgridppt) sp_count_i
550         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
551         write(unitoutgridppt) sp_count_r
552         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
553!         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
554
555
556  ! Mixing ratios
557
558  ! surf_only write only 1st layer
559
560         sp_count_i=0
561         sp_count_r=0
562         sp_fact=-1.
563         sp_zer=.true.
564          do kz=1,1
565            do jy=0,numygrid-1
566              do ix=0,numxgrid-1
567                if (grid(ix,jy,kz).gt.smallnum) then
568                  if (sp_zer.eqv..true.) then ! first non zero value
569                    sp_count_i=sp_count_i+1
570                    sparse_dump_i(sp_count_i)= &
571                         ix+jy*numxgrid+kz*numxgrid*numygrid
572                    sp_zer=.false.
573                    sp_fact=sp_fact*(-1.)
574                 endif
575                 sp_count_r=sp_count_r+1
576                 sparse_dump_r(sp_count_r)= &
577                      sp_fact* &
578                      1.e12*grid(ix,jy,kz) &
579                      /volume(ix,jy,kz)/outnum* &
580                      weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
581                 sparse_dump_u(sp_count_r)= &
582                      1.e12*gridsigma(ix,jy,kz)/volume(ix,jy,kz)/ &
583                      outnum*weightair/weightmolar(ks)/ &
584                      densityoutgrid(ix,jy,kz)
585              else ! concentration is zero
586                  sp_zer=.true.
587              endif
588              end do
589            end do
590          end do
591         write(unitoutgridppt) sp_count_i
592         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
593         write(unitoutgridppt) sp_count_r
594         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
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_surf
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG