source: flexpart.git/src/concoutput_surf_nest_mpi.f90 @ 4c64400

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

Parallelization of domain fill option (save/restart not implemented yet)

  • Property mode set to 100644
File size: 24.7 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_nest(itime,outnum)
23  !                        i     i
24  !*****************************************************************************
25  !                                                                            *
26  !     Output of the concentration grid and the receptor concentrations.      *
27  !                                                                            *
28  !     Author: A. Stohl                                                       *
29  !                                                                            *
30  !     24 May 1995                                                            *
31  !                                                                            *
32  !     13 April 1999, Major update: if output size is smaller, dump output    *
33  !                    in sparse matrix format; additional output of           *
34  !                    uncertainty                                             *
35  !                                                                            *
36  !     05 April 2000, Major update: output of age classes; output for backward*
37  !                    runs is time spent in grid cell times total mass of     *
38  !                    species.                                                *
39  !                                                                            *
40  !     17 February 2002, Appropriate dimensions for backward and forward runs *
41  !                       are now specified in file par_mod                    *
42  !                                                                            *
43  !     June 2006, write grid in sparse matrix with a single write command     *
44  !                in order to save disk space                                 *
45  !                                                                            *
46  !     2008 new sparse matrix format                                          *
47  !                                                                            *
48  !     2014 eso: MPI version. Only called by root process                     *
49  !                                                                            *
50  !*****************************************************************************
51  !                                                                            *
52  ! Variables:                                                                 *
53  ! outnum          number of samples                                          *
54  ! ncells          number of cells with non-zero concentrations               *
55  ! sparse          .true. if in sparse matrix format, else .false.            *
56  ! tot_mu          1 for forward, initial mass mixing ration for backw. runs  *
57  !                                                                            *
58  !*****************************************************************************
59
60  use unc_mod
61  use point_mod
62  use outg_mod
63  use par_mod
64  use com_mod
65  use mpi_mod
66  use mean_mod
67
68  implicit none
69
70  real(kind=dp) :: jul
71  integer :: itime,i,ix,jy,kz,ks,kp,l,iix,jjy,kzz,nage,jjjjmmdd,ihmmss
72  integer :: sp_count_i,sp_count_r
73  real :: sp_fact
74  real :: outnum,densityoutrecept(maxreceptor),xl,yl
75
76  !real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
77  !    +grid(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,maxpointspec_act,
78  !    +    maxageclass)
79  !real wetgrid(0:numxgrid-1,0:numygrid-1,maxspec,maxpointspec_act,
80  !    +       maxageclass)
81  !real drygrid(0:numxgrid-1,0:numygrid-1,maxspec,
82  !    +       maxpointspec_act,maxageclass)
83  !real gridsigma(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,
84  !    +       maxpointspec_act,maxageclass),
85  !    +     drygridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
86  !    +     maxpointspec_act,maxageclass),
87  !    +     wetgridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
88  !    +     maxpointspec_act,maxageclass)
89  !real factor(0:numxgrid-1,0:numygrid-1,numzgrid)
90  !real sparse_dump_r(numxgrid*numygrid*numzgrid)
91  !integer sparse_dump_i(numxgrid*numygrid*numzgrid)
92
93  !real sparse_dump_u(numxgrid*numygrid*numzgrid)
94  real(dep_prec) :: auxgrid(nclassunc)
95  real :: halfheight,dz,dz1,dz2,tot_mu(maxspec,maxpointspec_act)
96  real,parameter :: smallnum = tiny(0.0) ! smallest number that can be handled
97  real,parameter :: weightair=28.97
98  logical :: sp_zer
99  character :: adate*8,atime*6
100  character(len=3) :: anspec
101  integer :: mind
102! mind        eso:added to get consistent results between 2&3-fields versions
103
104! Measure execution time
105  if (mp_measure_time) call mpif_mtime('iotime',0)
106  !   call cpu_time(mp_root_time_beg)
107  !   mp_root_wtime_beg = mpi_wtime()
108  ! end if
109
110  if (verbosity.eq.1) then
111     print*,'inside concoutput_surf '
112     CALL SYSTEM_CLOCK(count_clock)
113     WRITE(*,*) 'SYSTEM_CLOCK',count_clock - count_clock0   
114  endif
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
124
125  ! For forward simulations, output fields have dimension MAXSPEC,
126  ! for backward simulations, output fields have dimension MAXPOINT.
127  ! Thus, make loops either about nspec, or about numpoint
128  !*****************************************************************
129
130
131    if (ldirect.eq.1) then
132       do ks=1,nspec
133         do kp=1,maxpointspec_act
134           tot_mu(ks,kp)=1
135         end do
136       end do
137   else
138      do ks=1,nspec
139             do kp=1,maxpointspec_act
140               tot_mu(ks,kp)=xmass(kp,ks)
141             end do
142      end do
143    endif
144
145
146  !*******************************************************************
147  ! Compute air density: sufficiently accurate to take it
148  ! from coarse grid at some time
149  ! Determine center altitude of output layer, and interpolate density
150  ! data to that altitude
151  !*******************************************************************
152
153  mind=memind(2)
154  do kz=1,numzgrid
155    if (kz.eq.1) then
156      halfheight=outheight(1)/2.
157    else
158      halfheight=(outheight(kz)+outheight(kz-1))/2.
159    endif
160    do kzz=2,nz
161      if ((height(kzz-1).lt.halfheight).and. &
162           (height(kzz).gt.halfheight)) goto 46
163    end do
16446   kzz=max(min(kzz,nz),2)
165    dz1=halfheight-height(kzz-1)
166    dz2=height(kzz)-halfheight
167    dz=dz1+dz2
168    do jy=0,numygridn-1
169      do ix=0,numxgridn-1
170        xl=outlon0n+real(ix)*dxoutn
171        yl=outlat0n+real(jy)*dyoutn
172        xl=(xl-xlon0)/dx
173        yl=(yl-ylat0)/dy
174        iix=max(min(nint(xl),nxmin1),0)
175        jjy=max(min(nint(yl),nymin1),0)
176        ! densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
177        !      rho(iix,jjy,kzz-1,2)*dz2)/dz
178        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,mind)*dz1+ &
179             rho(iix,jjy,kzz-1,mind)*dz2)/dz
180      end do
181    end do
182  end do
183
184  do i=1,numreceptor
185    xl=xreceptor(i)
186    yl=yreceptor(i)
187    iix=max(min(nint(xl),nxmin1),0)
188    jjy=max(min(nint(yl),nymin1),0)
189    !densityoutrecept(i)=rho(iix,jjy,1,2)
190    densityoutrecept(i)=rho(iix,jjy,1,mind)
191  end do
192
193
194  ! Output is different for forward and backward simulations
195    do kz=1,numzgrid
196      do jy=0,numygridn-1
197        do ix=0,numxgridn-1
198          if (ldirect.eq.1) then
199            factor3d(ix,jy,kz)=1.e12/volumen(ix,jy,kz)/outnum
200          else
201            factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
202          endif
203        end do
204      end do
205    end do
206
207  !*********************************************************************
208  ! Determine the standard deviation of the mean concentration or mixing
209  ! ratio (uncertainty of the output) and the dry and wet deposition
210  !*********************************************************************
211
212  do ks=1,nspec
213
214  write(anspec,'(i3.3)') ks
215  if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
216    if (ldirect.eq.1) then
217      open(unitoutgrid,file=path(2)(1:length(2))//'grid_conc_nest_' &
218           //adate// &
219           atime//'_'//anspec,form='unformatted')
220    else
221      open(unitoutgrid,file=path(2)(1:length(2))//'grid_time_nest_' &
222           //adate// &
223           atime//'_'//anspec,form='unformatted')
224    endif
225     write(unitoutgrid) itime
226   endif
227
228  if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
229   open(unitoutgridppt,file=path(2)(1:length(2))//'grid_pptv_nest_' &
230        //adate// &
231        atime//'_'//anspec,form='unformatted')
232
233    write(unitoutgridppt) itime
234  endif
235
236  do kp=1,maxpointspec_act
237  do nage=1,nageclass
238
239    do jy=0,numygridn-1
240      do ix=0,numxgridn-1
241
242  ! WET DEPOSITION
243        if ((WETDEP).and.(ldirect.gt.0)) then
244            do l=1,nclassunc
245              auxgrid(l)=wetgriduncn0(ix,jy,ks,kp,l,nage)
246            end do
247            call mean(auxgrid,wetgrid(ix,jy), &
248                 wetgridsigma(ix,jy),nclassunc)
249  ! Multiply by number of classes to get total concentration
250            wetgrid(ix,jy)=wetgrid(ix,jy) &
251                 *nclassunc
252  ! Calculate standard deviation of the mean
253            wetgridsigma(ix,jy)= &
254                 wetgridsigma(ix,jy)* &
255                 sqrt(real(nclassunc))
256        endif
257
258  ! DRY DEPOSITION
259        if ((DRYDEP).and.(ldirect.gt.0)) then
260            do l=1,nclassunc
261              auxgrid(l)=drygriduncn0(ix,jy,ks,kp,l,nage)
262            end do
263            call mean(auxgrid,drygrid(ix,jy), &
264                 drygridsigma(ix,jy),nclassunc)
265  ! Multiply by number of classes to get total concentration
266            drygrid(ix,jy)=drygrid(ix,jy)* &
267                 nclassunc
268  ! Calculate standard deviation of the mean
269            drygridsigma(ix,jy)= &
270                 drygridsigma(ix,jy)* &
271                 sqrt(real(nclassunc))
272        endif
273
274  ! CONCENTRATION OR MIXING RATIO
275        do kz=1,numzgrid
276            do l=1,nclassunc
277              auxgrid(l)=griduncn(ix,jy,kz,ks,kp,l,nage)
278            end do
279            call mean(auxgrid,grid(ix,jy,kz), &
280                 gridsigma(ix,jy,kz),nclassunc)
281  ! Multiply by number of classes to get total concentration
282            grid(ix,jy,kz)= &
283                 grid(ix,jy,kz)*nclassunc
284  ! Calculate standard deviation of the mean
285            gridsigma(ix,jy,kz)= &
286                 gridsigma(ix,jy,kz)* &
287                 sqrt(real(nclassunc))
288        end do
289      end do
290    end do
291
292
293  !*******************************************************************
294  ! Generate output: may be in concentration (ng/m3) or in mixing
295  ! ratio (ppt) or both
296  ! Output the position and the values alternated multiplied by
297  ! 1 or -1, first line is number of values, number of positions
298  ! For backward simulations, the unit is seconds, stored in grid_time
299  !*******************************************************************
300
301  ! Concentration output
302  !*********************
303  if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
304
305  ! Wet deposition
306         sp_count_i=0
307         sp_count_r=0
308         sp_fact=-1.
309         sp_zer=.true.
310         if ((ldirect.eq.1).and.(WETDEP)) then
311         do jy=0,numygridn-1
312            do ix=0,numxgridn-1
313  !oncentraion greater zero
314              if (wetgrid(ix,jy).gt.smallnum) then
315                 if (sp_zer.eqv..true.) then ! first non zero value
316                    sp_count_i=sp_count_i+1
317                    sparse_dump_i(sp_count_i)=ix+jy*numxgridn
318                    sp_zer=.false.
319                    sp_fact=sp_fact*(-1.)
320                 endif
321                 sp_count_r=sp_count_r+1
322                 sparse_dump_r(sp_count_r)= &
323                      sp_fact*1.e12*wetgrid(ix,jy)/arean(ix,jy)
324                 sparse_dump_u(sp_count_r)= &
325                      1.e12*wetgridsigma(ix,jy)/area(ix,jy)
326              else ! concentration is zero
327                  sp_zer=.true.
328              endif
329            end do
330         end do
331         else
332            sp_count_i=0
333            sp_count_r=0
334         endif
335         write(unitoutgrid) sp_count_i
336         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
337         write(unitoutgrid) sp_count_r
338         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
339         write(unitoutgrid) sp_count_r
340         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
341
342  ! Dry deposition
343         sp_count_i=0
344         sp_count_r=0
345         sp_fact=-1.
346         sp_zer=.true.
347         if ((ldirect.eq.1).and.(DRYDEP)) then
348          do jy=0,numygridn-1
349            do ix=0,numxgridn-1
350              if (drygrid(ix,jy).gt.smallnum) then
351                 if (sp_zer.eqv..true.) then ! first non zero value
352                    sp_count_i=sp_count_i+1
353                    sparse_dump_i(sp_count_i)=ix+jy*numxgridn
354                    sp_zer=.false.
355                    sp_fact=sp_fact*(-1.)
356                 endif
357                 sp_count_r=sp_count_r+1
358                 sparse_dump_r(sp_count_r)= &
359                      sp_fact* &
360                      1.e12*drygrid(ix,jy)/arean(ix,jy)
361                 sparse_dump_u(sp_count_r)= &
362                      1.e12*drygridsigma(ix,jy)/area(ix,jy)
363              else ! concentration is zero
364                  sp_zer=.true.
365              endif
366            end do
367          end do
368         else
369            sp_count_i=0
370            sp_count_r=0
371         endif
372         write(unitoutgrid) sp_count_i
373         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
374         write(unitoutgrid) sp_count_r
375         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
376         write(unitoutgrid) sp_count_r
377         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
378
379
380
381  ! Concentrations
382
383  ! if surf_only write only 1st layer
384
385         if(surf_only.eq.1) then
386         sp_count_i=0
387         sp_count_r=0
388         sp_fact=-1.
389         sp_zer=.true.
390          do kz=1,1
391            do jy=0,numygridn-1
392              do ix=0,numxgridn-1
393                if (grid(ix,jy,kz).gt.smallnum) then
394                  if (sp_zer.eqv..true.) then ! first non zero value
395                    sp_count_i=sp_count_i+1
396                    sparse_dump_i(sp_count_i)= &
397                         ix+jy*numxgridn+kz*numxgridn*numygridn
398                    sp_zer=.false.
399                    sp_fact=sp_fact*(-1.)
400                   endif
401                   sp_count_r=sp_count_r+1
402                   sparse_dump_r(sp_count_r)= &
403                        sp_fact* &
404                        grid(ix,jy,kz)* &
405                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
406  !                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
407  !    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
408                   sparse_dump_u(sp_count_r)= &
409                        gridsigma(ix,jy,kz)* &
410                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
411              else ! concentration is zero
412                  sp_zer=.true.
413              endif
414              end do
415            end do
416          end do
417         write(unitoutgrid) sp_count_i
418         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
419         write(unitoutgrid) sp_count_r
420         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
421         write(unitoutgrid) sp_count_r
422         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
423         else
424
425  ! write full vertical resolution
426
427         sp_count_i=0
428         sp_count_r=0
429         sp_fact=-1.
430         sp_zer=.true.
431          do kz=1,numzgrid
432            do jy=0,numygridn-1
433              do ix=0,numxgridn-1
434                if (grid(ix,jy,kz).gt.smallnum) then
435                  if (sp_zer.eqv..true.) then ! first non zero value
436                    sp_count_i=sp_count_i+1
437                    sparse_dump_i(sp_count_i)= &
438                         ix+jy*numxgridn+kz*numxgridn*numygridn
439                    sp_zer=.false.
440                    sp_fact=sp_fact*(-1.)
441                   endif
442                   sp_count_r=sp_count_r+1
443                   sparse_dump_r(sp_count_r)= &
444                        sp_fact* &
445                        grid(ix,jy,kz)* &
446                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
447  !                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
448  !    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
449                   sparse_dump_u(sp_count_r)= &
450                        gridsigma(ix,jy,kz)* &
451                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
452              else ! concentration is zero
453                  sp_zer=.true.
454              endif
455              end do
456            end do
457          end do
458         write(unitoutgrid) sp_count_i
459         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
460         write(unitoutgrid) sp_count_r
461         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
462         write(unitoutgrid) sp_count_r
463         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
464         endif ! surf_only
465
466
467    endif !  concentration output
468
469  ! Mixing ratio output
470  !********************
471
472  if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
473
474  ! Wet deposition
475         sp_count_i=0
476         sp_count_r=0
477         sp_fact=-1.
478         sp_zer=.true.
479         if ((ldirect.eq.1).and.(WETDEP)) then
480          do jy=0,numygridn-1
481            do ix=0,numxgridn-1
482                if (wetgrid(ix,jy).gt.smallnum) then
483                  if (sp_zer.eqv..true.) then ! first non zero value
484                    sp_count_i=sp_count_i+1
485                    sparse_dump_i(sp_count_i)= &
486                         ix+jy*numxgridn
487                    sp_zer=.false.
488                    sp_fact=sp_fact*(-1.)
489                 endif
490                 sp_count_r=sp_count_r+1
491                 sparse_dump_r(sp_count_r)= &
492                      sp_fact* &
493                      1.e12*wetgrid(ix,jy)/arean(ix,jy)
494                 sparse_dump_u(sp_count_r)= &
495                      1.e12*wetgridsigma(ix,jy)/area(ix,jy)
496              else ! concentration is zero
497                  sp_zer=.true.
498              endif
499            end do
500          end do
501         else
502           sp_count_i=0
503           sp_count_r=0
504         endif
505         write(unitoutgridppt) sp_count_i
506         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
507         write(unitoutgridppt) sp_count_r
508         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
509         write(unitoutgridppt) sp_count_r
510         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
511
512
513  ! Dry deposition
514         sp_count_i=0
515         sp_count_r=0
516         sp_fact=-1.
517         sp_zer=.true.
518         if ((ldirect.eq.1).and.(DRYDEP)) then
519          do jy=0,numygridn-1
520            do ix=0,numxgridn-1
521                if (drygrid(ix,jy).gt.smallnum) then
522                  if (sp_zer.eqv..true.) then ! first non zero value
523                    sp_count_i=sp_count_i+1
524                    sparse_dump_i(sp_count_i)= &
525                         ix+jy*numxgridn
526                    sp_zer=.false.
527                    sp_fact=sp_fact*(-1)
528                 endif
529                 sp_count_r=sp_count_r+1
530                 sparse_dump_r(sp_count_r)= &
531                      sp_fact* &
532                      1.e12*drygrid(ix,jy)/arean(ix,jy)
533                 sparse_dump_u(sp_count_r)= &
534                      1.e12*drygridsigma(ix,jy)/area(ix,jy)
535              else ! concentration is zero
536                  sp_zer=.true.
537              endif
538            end do
539          end do
540         else
541           sp_count_i=0
542           sp_count_r=0
543         endif
544         write(unitoutgridppt) sp_count_i
545         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
546         write(unitoutgridppt) sp_count_r
547         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
548         write(unitoutgridppt) sp_count_r
549         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
550
551
552  ! Mixing ratios
553
554    ! if surf_only write only 1st layer
555
556         if(surf_only.eq.1) then
557         sp_count_i=0
558         sp_count_r=0
559         sp_fact=-1.
560         sp_zer=.true.
561          do kz=1,1
562            do jy=0,numygridn-1
563              do ix=0,numxgridn-1
564                if (grid(ix,jy,kz).gt.smallnum) then
565                  if (sp_zer.eqv..true.) then ! first non zero value
566                    sp_count_i=sp_count_i+1
567                    sparse_dump_i(sp_count_i)= &
568                         ix+jy*numxgridn+kz*numxgridn*numygridn
569                    sp_zer=.false.
570                    sp_fact=sp_fact*(-1.)
571                 endif
572                 sp_count_r=sp_count_r+1
573                 sparse_dump_r(sp_count_r)= &
574                      sp_fact* &
575                      1.e12*grid(ix,jy,kz) &
576                      /volumen(ix,jy,kz)/outnum* &
577                      weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
578                 sparse_dump_u(sp_count_r)= &
579                      1.e12*gridsigma(ix,jy,kz)/volumen(ix,jy,kz)/ &
580                      outnum*weightair/weightmolar(ks)/ &
581                      densityoutgrid(ix,jy,kz)
582              else ! concentration is zero
583                  sp_zer=.true.
584              endif
585              end do
586            end do
587          end do
588         write(unitoutgridppt) sp_count_i
589         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
590         write(unitoutgridppt) sp_count_r
591         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
592         write(unitoutgridppt) sp_count_r
593         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
594         else
595
596  ! write full vertical resolution
597
598         sp_count_i=0
599         sp_count_r=0
600         sp_fact=-1.
601         sp_zer=.true.
602          do kz=1,numzgrid
603            do jy=0,numygridn-1
604              do ix=0,numxgridn-1
605                if (grid(ix,jy,kz).gt.smallnum) then
606                  if (sp_zer.eqv..true.) then ! first non zero value
607                    sp_count_i=sp_count_i+1
608                    sparse_dump_i(sp_count_i)= &
609                         ix+jy*numxgridn+kz*numxgridn*numygridn
610                    sp_zer=.false.
611                    sp_fact=sp_fact*(-1.)
612                 endif
613                 sp_count_r=sp_count_r+1
614                 sparse_dump_r(sp_count_r)= &
615                      sp_fact* &
616                      1.e12*grid(ix,jy,kz) &
617                      /volumen(ix,jy,kz)/outnum* &
618                      weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
619                 sparse_dump_u(sp_count_r)= &
620                      1.e12*gridsigma(ix,jy,kz)/volumen(ix,jy,kz)/ &
621                      outnum*weightair/weightmolar(ks)/ &
622                      densityoutgrid(ix,jy,kz)
623              else ! concentration is zero
624                  sp_zer=.true.
625              endif
626              end do
627            end do
628          end do
629         write(unitoutgridppt) sp_count_i
630         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
631         write(unitoutgridppt) sp_count_r
632         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
633         write(unitoutgridppt) sp_count_r
634         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
635         endif ! surf_only
636
637      endif ! output for ppt
638
639  end do
640  end do
641
642    close(unitoutgridppt)
643    close(unitoutgrid)
644
645  end do
646
647
648
649  ! Reinitialization of grid
650  !*************************
651
652  do ks=1,nspec
653  do kp=1,maxpointspec_act
654    do i=1,numreceptor
655      creceptor(i,ks)=0.
656    end do
657    do jy=0,numygridn-1
658      do ix=0,numxgridn-1
659        do l=1,nclassunc
660          do nage=1,nageclass
661            do kz=1,numzgrid
662              griduncn(ix,jy,kz,ks,kp,l,nage)=0.
663            end do
664          end do
665        end do
666      end do
667    end do
668  end do
669  end do
670
671  if (mp_measure_time) call mpif_mtime('iotime',1)
672  ! if (mp_measure_time) then
673  !   call cpu_time(mp_root_time_end)
674  !   mp_root_wtime_end = mpi_wtime()
675  !   mp_root_time_total = mp_root_time_total + (mp_root_time_end - mp_root_time_beg)
676  !   mp_root_wtime_total = mp_root_wtime_total + (mp_root_wtime_end - mp_root_wtime_beg)
677  ! end if
678
679end subroutine concoutput_surf_nest
680
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG