Changeset b069789 in flexpart.git


Ignore:
Timestamp:
Oct 16, 2015, 11:31:33 AM (8 years ago)
Author:
Espen Sollum ATMOS <eso@…>
Branches:
master, 10.4.1_pesei, GFS_025, bugfixes+enhancements, dev, release-10, release-10.4.1, scaling-bug, univie
Children:
c04b739
Parents:
43225d1
Message:

If the 'dates' file exists on simulation start, backup to 'dates.bak' and create new file 'dates' instead of appending to it

Location:
src
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • src/concoutput.f90

    r478e9e6 rb069789  
    2121
    2222subroutine 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   !*****************************************************************************
    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   !*****************************************************************************
     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!*****************************************************************************
     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!*****************************************************************************
    5959
    6060  use unc_mod
     
    7272  real :: outnum,densityoutrecept(maxreceptor),xl,yl
    7373
    74   !real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
    75   !    +grid(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,maxpointspec_act,
    76   !    +    maxageclass)
    77   !real wetgrid(0:numxgrid-1,0:numygrid-1,maxspec,maxpointspec_act,
    78   !    +       maxageclass)
    79   !real drygrid(0:numxgrid-1,0:numygrid-1,maxspec,
    80   !    +       maxpointspec_act,maxageclass)
    81   !real gridsigma(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,
    82   !    +       maxpointspec_act,maxageclass),
    83   !    +     drygridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
    84   !    +     maxpointspec_act,maxageclass),
    85   !    +     wetgridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
    86   !    +     maxpointspec_act,maxageclass)
    87   !real factor(0:numxgrid-1,0:numygrid-1,numzgrid)
    88   !real sparse_dump_r(numxgrid*numygrid*numzgrid)
    89   !integer sparse_dump_i(numxgrid*numygrid*numzgrid)
    90 
    91   !real sparse_dump_u(numxgrid*numygrid*numzgrid)
     74!real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
     75!    +grid(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,maxpointspec_act,
     76!    +    maxageclass)
     77!real wetgrid(0:numxgrid-1,0:numygrid-1,maxspec,maxpointspec_act,
     78!    +       maxageclass)
     79!real drygrid(0:numxgrid-1,0:numygrid-1,maxspec,
     80!    +       maxpointspec_act,maxageclass)
     81!real gridsigma(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,
     82!    +       maxpointspec_act,maxageclass),
     83!    +     drygridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
     84!    +     maxpointspec_act,maxageclass),
     85!    +     wetgridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
     86!    +     maxpointspec_act,maxageclass)
     87!real factor(0:numxgrid-1,0:numygrid-1,numzgrid)
     88!real sparse_dump_r(numxgrid*numygrid*numzgrid)
     89!integer sparse_dump_i(numxgrid*numygrid*numzgrid)
     90
     91!real sparse_dump_u(numxgrid*numygrid*numzgrid)
    9292  real :: auxgrid(nclassunc),gridtotal,gridsigmatotal,gridtotalunc
    9393  real :: wetgridtotal,wetgridsigmatotal,wetgridtotalunc
     
    9797  real,parameter :: weightair=28.97
    9898  logical :: sp_zer
    99   LOGICAL,save :: init=.true.
     99  logical,save :: init=.true.
    100100  character :: adate*8,atime*6
    101101  character(len=3) :: anspec
    102102  integer :: mind
    103103! mind        eso:added to ensure identical results between 2&3-fields versions
    104   CHARACTER(LEN=8),save :: file_stat='REPLACE'
    105 
    106   ! Determine current calendar date, needed for the file name
    107   !**********************************************************
     104  character(LEN=8),save :: file_stat='REPLACE'
     105  logical :: ldates_file
     106  integer :: ierr
     107  character(LEN=100) :: dates_char
     108
     109! Determine current calendar date, needed for the file name
     110!**********************************************************
    108111
    109112  jul=bdate+real(itime,kind=dp)/86400._dp
     
    111114  write(adate,'(i8.8)') jjjjmmdd
    112115  write(atime,'(i6.6)') ihmmss
    113   OPEN(unitdates,file=path(2)(1:length(2))//'dates', ACCESS='APPEND', STATUS=file_stat)
     116
     117! Overwrite existing dates file on first call, later append to it
     118! This fixes a bug where the dates file kept growing across multiple runs
     119
     120! If 'dates' file exists, make a backup
     121  inquire(file=path(2)(1:length(2))//'dates', exist=ldates_file)
     122  if (ldates_file.and.init) then
     123    open(unit=unitdates, file=path(2)(1:length(2))//'dates',form='formatted', &
     124         &access='sequential', status='old', action='read', iostat=ierr)
     125    open(unit=unittmp, file=path(2)(1:length(2))//'dates.bak', access='sequential', &
     126         &status='replace', action='write', form='formatted', iostat=ierr)
     127    do while (.true.)
     128      read(unitdates, '(a)', iostat=ierr) dates_char
     129      if (ierr.ne.0) exit
     130      write(unit=unittmp, fmt='(a)', iostat=ierr, advance='yes') trim(dates_char)
     131!      if (ierr.ne.0) write(*,*) "Write error, ", ierr
     132    end do
     133    close(unit=unitdates)
     134    close(unit=unittmp)
     135  end if
     136
     137  open(unitdates,file=path(2)(1:length(2))//'dates', ACCESS='APPEND', STATUS=file_stat)
    114138  write(unitdates,'(a)') adate//atime
    115139  close(unitdates) 
    116140
    117   ! Overwrite existing dates file on first call, later append to it
    118   ! This fixes a bug where the dates file kept growing across multiple runs
    119   ! TODO check if the 'always APPEND'-behaviour is useful in other scenarioes
    120   ! e.g. (restart?)
    121141  IF (init) THEN
    122142    file_stat='OLD'
     
    125145
    126146
    127 
    128   ! For forward simulations, output fields have dimension MAXSPEC,
    129   ! for backward simulations, output fields have dimension MAXPOINT.
    130   ! Thus, make loops either about nspec, or about numpoint
    131   !*****************************************************************
     147! For forward simulations, output fields have dimension MAXSPEC,
     148! for backward simulations, output fields have dimension MAXPOINT.
     149! Thus, make loops either about nspec, or about numpoint
     150!*****************************************************************
    132151
    133152
     
    147166
    148167
    149   !*******************************************************************
    150   ! Compute air density: sufficiently accurate to take it
    151   ! from coarse grid at some time
    152   ! Determine center altitude of output layer, and interpolate density
    153   ! data to that altitude
    154   !*******************************************************************
     168!*******************************************************************
     169! Compute air density: sufficiently accurate to take it
     170! from coarse grid at some time
     171! Determine center altitude of output layer, and interpolate density
     172! data to that altitude
     173!*******************************************************************
    155174
    156175  mind=memind(2)
     
    165184           (height(kzz).gt.halfheight)) goto 46
    166185    end do
    167 46   kzz=max(min(kzz,nz),2)
     18646  kzz=max(min(kzz,nz),2)
    168187    dz1=halfheight-height(kzz-1)
    169188    dz2=height(kzz)-halfheight
     
    177196        iix=max(min(nint(xl),nxmin1),0)
    178197        jjy=max(min(nint(yl),nymin1),0)
    179         ! densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
    180         !      rho(iix,jjy,kzz-1,2)*dz2)/dz
     198! densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
     199!      rho(iix,jjy,kzz-1,2)*dz2)/dz
    181200        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,mind)*dz1+ &
    182201             rho(iix,jjy,kzz-1,mind)*dz2)/dz
     
    190209    iix=max(min(nint(xl),nxmin1),0)
    191210    jjy=max(min(nint(yl),nymin1),0)
    192     !densityoutrecept(i)=rho(iix,jjy,1,2)
     211!densityoutrecept(i)=rho(iix,jjy,1,2)
    193212    densityoutrecept(i)=rho(iix,jjy,1,mind)
    194213  end do
    195214
    196215
    197   ! Output is different for forward and backward simulations
    198     do kz=1,numzgrid
    199       do jy=0,numygrid-1
    200         do ix=0,numxgrid-1
    201           if (ldirect.eq.1) then
    202             factor3d(ix,jy,kz)=1.e12/volume(ix,jy,kz)/outnum
    203           else
    204             factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
    205           endif
    206         end do
     216! Output is different for forward and backward simulations
     217  do kz=1,numzgrid
     218    do jy=0,numygrid-1
     219      do ix=0,numxgrid-1
     220        if (ldirect.eq.1) then
     221          factor3d(ix,jy,kz)=1.e12/volume(ix,jy,kz)/outnum
     222        else
     223          factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
     224        endif
    207225      end do
    208226    end do
    209 
    210   !*********************************************************************
    211   ! Determine the standard deviation of the mean concentration or mixing
    212   ! ratio (uncertainty of the output) and the dry and wet deposition
    213   !*********************************************************************
     227  end do
     228
     229!*********************************************************************
     230! Determine the standard deviation of the mean concentration or mixing
     231! ratio (uncertainty of the output) and the dry and wet deposition
     232!*********************************************************************
    214233
    215234  gridtotal=0.
     
    225244  do ks=1,nspec
    226245
    227   write(anspec,'(i3.3)') ks
    228   if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
    229     if (ldirect.eq.1) then
    230       open(unitoutgrid,file=path(2)(1:length(2))//'grid_conc_'//adate// &
     246    write(anspec,'(i3.3)') ks
     247    if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
     248      if (ldirect.eq.1) then
     249        open(unitoutgrid,file=path(2)(1:length(2))//'grid_conc_'//adate// &
     250             atime//'_'//anspec,form='unformatted')
     251      else
     252        open(unitoutgrid,file=path(2)(1:length(2))//'grid_time_'//adate// &
     253             atime//'_'//anspec,form='unformatted')
     254      endif
     255      write(unitoutgrid) itime
     256    endif
     257
     258    if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
     259      open(unitoutgridppt,file=path(2)(1:length(2))//'grid_pptv_'//adate// &
    231260           atime//'_'//anspec,form='unformatted')
    232     else
    233       open(unitoutgrid,file=path(2)(1:length(2))//'grid_time_'//adate// &
    234            atime//'_'//anspec,form='unformatted')
     261
     262      write(unitoutgridppt) itime
    235263    endif
    236      write(unitoutgrid) itime
    237    endif
    238 
    239   if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
    240    open(unitoutgridppt,file=path(2)(1:length(2))//'grid_pptv_'//adate// &
    241         atime//'_'//anspec,form='unformatted')
    242 
    243     write(unitoutgridppt) itime
    244   endif
    245 
    246   do kp=1,maxpointspec_act
    247   do nage=1,nageclass
    248 
    249     do jy=0,numygrid-1
    250       do ix=0,numxgrid-1
    251 
    252   ! WET DEPOSITION
    253         if ((WETDEP).and.(ldirect.gt.0)) then
    254             do l=1,nclassunc
    255               auxgrid(l)=wetgridunc(ix,jy,ks,kp,l,nage)
    256             end do
    257             call mean(auxgrid,wetgrid(ix,jy), &
    258                  wetgridsigma(ix,jy),nclassunc)
    259   ! Multiply by number of classes to get total concentration
    260             wetgrid(ix,jy)=wetgrid(ix,jy) &
    261                  *nclassunc
    262             wetgridtotal=wetgridtotal+wetgrid(ix,jy)
    263   ! Calculate standard deviation of the mean
    264             wetgridsigma(ix,jy)= &
    265                  wetgridsigma(ix,jy)* &
    266                  sqrt(real(nclassunc))
    267             wetgridsigmatotal=wetgridsigmatotal+ &
    268                  wetgridsigma(ix,jy)
    269         endif
    270 
    271   ! DRY DEPOSITION
    272         if ((DRYDEP).and.(ldirect.gt.0)) then
    273             do l=1,nclassunc
    274               auxgrid(l)=drygridunc(ix,jy,ks,kp,l,nage)
    275             end do
    276             call mean(auxgrid,drygrid(ix,jy), &
    277                  drygridsigma(ix,jy),nclassunc)
    278   ! Multiply by number of classes to get total concentration
    279             drygrid(ix,jy)=drygrid(ix,jy)* &
    280                  nclassunc
    281             drygridtotal=drygridtotal+drygrid(ix,jy)
    282   ! Calculate standard deviation of the mean
    283             drygridsigma(ix,jy)= &
    284                  drygridsigma(ix,jy)* &
    285                  sqrt(real(nclassunc))
    286 125         drygridsigmatotal=drygridsigmatotal+ &
    287                  drygridsigma(ix,jy)
    288         endif
    289 
    290   ! CONCENTRATION OR MIXING RATIO
    291         do kz=1,numzgrid
    292             do l=1,nclassunc
    293               auxgrid(l)=gridunc(ix,jy,kz,ks,kp,l,nage)
    294             end do
    295             call mean(auxgrid,grid(ix,jy,kz), &
    296                  gridsigma(ix,jy,kz),nclassunc)
    297   ! Multiply by number of classes to get total concentration
    298             grid(ix,jy,kz)= &
    299                  grid(ix,jy,kz)*nclassunc
    300             gridtotal=gridtotal+grid(ix,jy,kz)
    301   ! Calculate standard deviation of the mean
    302             gridsigma(ix,jy,kz)= &
    303                  gridsigma(ix,jy,kz)* &
    304                  sqrt(real(nclassunc))
    305             gridsigmatotal=gridsigmatotal+ &
    306                  gridsigma(ix,jy,kz)
     264
     265    do kp=1,maxpointspec_act
     266      do nage=1,nageclass
     267
     268        do jy=0,numygrid-1
     269          do ix=0,numxgrid-1
     270
     271! WET DEPOSITION
     272            if ((WETDEP).and.(ldirect.gt.0)) then
     273              do l=1,nclassunc
     274                auxgrid(l)=wetgridunc(ix,jy,ks,kp,l,nage)
     275              end do
     276              call mean(auxgrid,wetgrid(ix,jy), &
     277                   wetgridsigma(ix,jy),nclassunc)
     278! Multiply by number of classes to get total concentration
     279              wetgrid(ix,jy)=wetgrid(ix,jy) &
     280                   *nclassunc
     281              wetgridtotal=wetgridtotal+wetgrid(ix,jy)
     282! Calculate standard deviation of the mean
     283              wetgridsigma(ix,jy)= &
     284                   wetgridsigma(ix,jy)* &
     285                   sqrt(real(nclassunc))
     286              wetgridsigmatotal=wetgridsigmatotal+ &
     287                   wetgridsigma(ix,jy)
     288            endif
     289
     290! DRY DEPOSITION
     291            if ((DRYDEP).and.(ldirect.gt.0)) then
     292              do l=1,nclassunc
     293                auxgrid(l)=drygridunc(ix,jy,ks,kp,l,nage)
     294              end do
     295              call mean(auxgrid,drygrid(ix,jy), &
     296                   drygridsigma(ix,jy),nclassunc)
     297! Multiply by number of classes to get total concentration
     298              drygrid(ix,jy)=drygrid(ix,jy)* &
     299                   nclassunc
     300              drygridtotal=drygridtotal+drygrid(ix,jy)
     301! Calculate standard deviation of the mean
     302              drygridsigma(ix,jy)= &
     303                   drygridsigma(ix,jy)* &
     304                   sqrt(real(nclassunc))
     305125           drygridsigmatotal=drygridsigmatotal+ &
     306                   drygridsigma(ix,jy)
     307            endif
     308
     309! CONCENTRATION OR MIXING RATIO
     310            do kz=1,numzgrid
     311              do l=1,nclassunc
     312                auxgrid(l)=gridunc(ix,jy,kz,ks,kp,l,nage)
     313              end do
     314              call mean(auxgrid,grid(ix,jy,kz), &
     315                   gridsigma(ix,jy,kz),nclassunc)
     316! Multiply by number of classes to get total concentration
     317              grid(ix,jy,kz)= &
     318                   grid(ix,jy,kz)*nclassunc
     319              gridtotal=gridtotal+grid(ix,jy,kz)
     320! Calculate standard deviation of the mean
     321              gridsigma(ix,jy,kz)= &
     322                   gridsigma(ix,jy,kz)* &
     323                   sqrt(real(nclassunc))
     324              gridsigmatotal=gridsigmatotal+ &
     325                   gridsigma(ix,jy,kz)
     326            end do
     327          end do
    307328        end do
    308       end do
    309     end do
    310 
    311 
    312 
    313 
    314   !*******************************************************************
    315   ! Generate output: may be in concentration (ng/m3) or in mixing
    316   ! ratio (ppt) or both
    317   ! Output the position and the values alternated multiplied by
    318   ! 1 or -1, first line is number of values, number of positions
    319   ! For backward simulations, the unit is seconds, stored in grid_time
    320   !*******************************************************************
    321 
    322   ! Concentration output
    323   !*********************
    324   if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
    325 
    326   ! Wet deposition
    327          sp_count_i=0
    328          sp_count_r=0
    329          sp_fact=-1.
    330          sp_zer=.true.
    331          if ((ldirect.eq.1).and.(WETDEP)) then
    332          do jy=0,numygrid-1
    333             do ix=0,numxgrid-1
    334   !oncentraion greater zero
    335               if (wetgrid(ix,jy).gt.smallnum) then
    336                  if (sp_zer.eqv..true.) then ! first non zero value
     329
     330
     331
     332
     333!*******************************************************************
     334! Generate output: may be in concentration (ng/m3) or in mixing
     335! ratio (ppt) or both
     336! Output the position and the values alternated multiplied by
     337! 1 or -1, first line is number of values, number of positions
     338! For backward simulations, the unit is seconds, stored in grid_time
     339!*******************************************************************
     340
     341! Concentration output
     342!*********************
     343        if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
     344
     345! Wet deposition
     346          sp_count_i=0
     347          sp_count_r=0
     348          sp_fact=-1.
     349          sp_zer=.true.
     350          if ((ldirect.eq.1).and.(WETDEP)) then
     351            do jy=0,numygrid-1
     352              do ix=0,numxgrid-1
     353!oncentraion greater zero
     354                if (wetgrid(ix,jy).gt.smallnum) then
     355                  if (sp_zer.eqv..true.) then ! first non zero value
    337356                    sp_count_i=sp_count_i+1
    338357                    sparse_dump_i(sp_count_i)=ix+jy*numxgrid
    339358                    sp_zer=.false.
    340359                    sp_fact=sp_fact*(-1.)
    341                  endif
    342                  sp_count_r=sp_count_r+1
    343                  sparse_dump_r(sp_count_r)= &
    344                       sp_fact*1.e12*wetgrid(ix,jy)/area(ix,jy)
    345   !                sparse_dump_u(sp_count_r)=
    346   !+                1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
    347               else ! concentration is zero
     360                  endif
     361                  sp_count_r=sp_count_r+1
     362                  sparse_dump_r(sp_count_r)= &
     363                       sp_fact*1.e12*wetgrid(ix,jy)/area(ix,jy)
     364!                sparse_dump_u(sp_count_r)=
     365!+                1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
     366                else ! concentration is zero
    348367                  sp_zer=.true.
    349               endif
    350             end do
    351          end do
    352          else
     368                endif
     369              end do
     370            end do
     371          else
    353372            sp_count_i=0
    354373            sp_count_r=0
    355          endif
    356          write(unitoutgrid) sp_count_i
    357          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
    358          write(unitoutgrid) sp_count_r
    359          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
    360   !       write(unitoutgrid) sp_count_u
    361   !       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
    362 
    363   ! Dry deposition
    364          sp_count_i=0
    365          sp_count_r=0
    366          sp_fact=-1.
    367          sp_zer=.true.
    368          if ((ldirect.eq.1).and.(DRYDEP)) then
    369           do jy=0,numygrid-1
    370             do ix=0,numxgrid-1
    371               if (drygrid(ix,jy).gt.smallnum) then
    372                  if (sp_zer.eqv..true.) then ! first non zero value
     374          endif
     375          write(unitoutgrid) sp_count_i
     376          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
     377          write(unitoutgrid) sp_count_r
     378          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
     379!       write(unitoutgrid) sp_count_u
     380!       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
     381
     382! Dry deposition
     383          sp_count_i=0
     384          sp_count_r=0
     385          sp_fact=-1.
     386          sp_zer=.true.
     387          if ((ldirect.eq.1).and.(DRYDEP)) then
     388            do jy=0,numygrid-1
     389              do ix=0,numxgrid-1
     390                if (drygrid(ix,jy).gt.smallnum) then
     391                  if (sp_zer.eqv..true.) then ! first non zero value
    373392                    sp_count_i=sp_count_i+1
    374393                    sparse_dump_i(sp_count_i)=ix+jy*numxgrid
    375394                    sp_zer=.false.
    376395                    sp_fact=sp_fact*(-1.)
    377                  endif
    378                  sp_count_r=sp_count_r+1
    379                  sparse_dump_r(sp_count_r)= &
    380                       sp_fact* &
    381                       1.e12*drygrid(ix,jy)/area(ix,jy)
    382   !                sparse_dump_u(sp_count_r)=
    383   !+                1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
    384               else ! concentration is zero
     396                  endif
     397                  sp_count_r=sp_count_r+1
     398                  sparse_dump_r(sp_count_r)= &
     399                       sp_fact* &
     400                       1.e12*drygrid(ix,jy)/area(ix,jy)
     401!                sparse_dump_u(sp_count_r)=
     402!+                1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
     403                else ! concentration is zero
    385404                  sp_zer=.true.
    386               endif
    387             end do
    388           end do
    389          else
     405                endif
     406              end do
     407            end do
     408          else
    390409            sp_count_i=0
    391410            sp_count_r=0
    392          endif
    393          write(unitoutgrid) sp_count_i
    394          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
    395          write(unitoutgrid) sp_count_r
    396          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
    397   !       write(*,*) sp_count_u
    398   !       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
    399 
    400 
    401 
    402   ! Concentrations
    403          sp_count_i=0
    404          sp_count_r=0
    405          sp_fact=-1.
    406          sp_zer=.true.
     411          endif
     412          write(unitoutgrid) sp_count_i
     413          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
     414          write(unitoutgrid) sp_count_r
     415          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
     416!       write(*,*) sp_count_u
     417!       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
     418
     419
     420
     421! Concentrations
     422          sp_count_i=0
     423          sp_count_r=0
     424          sp_fact=-1.
     425          sp_zer=.true.
    407426          do kz=1,numzgrid
    408427            do jy=0,numygrid-1
     
    415434                    sp_zer=.false.
    416435                    sp_fact=sp_fact*(-1.)
    417                    endif
    418                    sp_count_r=sp_count_r+1
    419                    sparse_dump_r(sp_count_r)= &
    420                         sp_fact* &
    421                         grid(ix,jy,kz)* &
    422                         factor3d(ix,jy,kz)/tot_mu(ks,kp)
    423   !                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
    424   !    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
    425   !                sparse_dump_u(sp_count_r)=
    426   !+               ,gridsigma(ix,jy,kz,ks,kp,nage)*
    427   !+               factor(ix,jy,kz)/tot_mu(ks,kp)
    428               else ! concentration is zero
     436                  endif
     437                  sp_count_r=sp_count_r+1
     438                  sparse_dump_r(sp_count_r)= &
     439                       sp_fact* &
     440                       grid(ix,jy,kz)* &
     441                       factor3d(ix,jy,kz)/tot_mu(ks,kp)
     442!                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
     443!    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
     444!                sparse_dump_u(sp_count_r)=
     445!+               ,gridsigma(ix,jy,kz,ks,kp,nage)*
     446!+               factor(ix,jy,kz)/tot_mu(ks,kp)
     447                else ! concentration is zero
    429448                  sp_zer=.true.
    430               endif
     449                endif
    431450              end do
    432451            end do
    433452          end do
    434          write(unitoutgrid) sp_count_i
    435          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
    436          write(unitoutgrid) sp_count_r
    437          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
    438   !       write(unitoutgrid) sp_count_u
    439   !       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
    440 
    441 
    442 
    443     endif !  concentration output
    444 
    445   ! Mixing ratio output
    446   !********************
    447 
    448   if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
    449 
    450   ! Wet deposition
    451          sp_count_i=0
    452          sp_count_r=0
    453          sp_fact=-1.
    454          sp_zer=.true.
    455          if ((ldirect.eq.1).and.(WETDEP)) then
    456           do jy=0,numygrid-1
    457             do ix=0,numxgrid-1
     453          write(unitoutgrid) sp_count_i
     454          write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
     455          write(unitoutgrid) sp_count_r
     456          write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
     457!       write(unitoutgrid) sp_count_u
     458!       write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
     459
     460
     461
     462        endif !  concentration output
     463
     464! Mixing ratio output
     465!********************
     466
     467        if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
     468
     469! Wet deposition
     470          sp_count_i=0
     471          sp_count_r=0
     472          sp_fact=-1.
     473          sp_zer=.true.
     474          if ((ldirect.eq.1).and.(WETDEP)) then
     475            do jy=0,numygrid-1
     476              do ix=0,numxgrid-1
    458477                if (wetgrid(ix,jy).gt.smallnum) then
    459478                  if (sp_zer.eqv..true.) then ! first non zero value
     
    463482                    sp_zer=.false.
    464483                    sp_fact=sp_fact*(-1.)
    465                  endif
    466                  sp_count_r=sp_count_r+1
    467                  sparse_dump_r(sp_count_r)= &
    468                       sp_fact* &
    469                       1.e12*wetgrid(ix,jy)/area(ix,jy)
    470   !                sparse_dump_u(sp_count_r)=
    471   !    +            ,1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
    472               else ! concentration is zero
     484                  endif
     485                  sp_count_r=sp_count_r+1
     486                  sparse_dump_r(sp_count_r)= &
     487                       sp_fact* &
     488                       1.e12*wetgrid(ix,jy)/area(ix,jy)
     489!                sparse_dump_u(sp_count_r)=
     490!    +            ,1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
     491                else ! concentration is zero
    473492                  sp_zer=.true.
    474               endif
    475             end do
    476           end do
    477          else
    478            sp_count_i=0
    479            sp_count_r=0
    480          endif
    481          write(unitoutgridppt) sp_count_i
    482          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
    483          write(unitoutgridppt) sp_count_r
    484          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
    485   !       write(unitoutgridppt) sp_count_u
    486   !       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
    487 
    488 
    489   ! Dry deposition
    490          sp_count_i=0
    491          sp_count_r=0
    492          sp_fact=-1.
    493          sp_zer=.true.
    494          if ((ldirect.eq.1).and.(DRYDEP)) then
    495           do jy=0,numygrid-1
    496             do ix=0,numxgrid-1
     493                endif
     494              end do
     495            end do
     496          else
     497            sp_count_i=0
     498            sp_count_r=0
     499          endif
     500          write(unitoutgridppt) sp_count_i
     501          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
     502          write(unitoutgridppt) sp_count_r
     503          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
     504!       write(unitoutgridppt) sp_count_u
     505!       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
     506
     507
     508! Dry deposition
     509          sp_count_i=0
     510          sp_count_r=0
     511          sp_fact=-1.
     512          sp_zer=.true.
     513          if ((ldirect.eq.1).and.(DRYDEP)) then
     514            do jy=0,numygrid-1
     515              do ix=0,numxgrid-1
    497516                if (drygrid(ix,jy).gt.smallnum) then
    498517                  if (sp_zer.eqv..true.) then ! first non zero value
     
    502521                    sp_zer=.false.
    503522                    sp_fact=sp_fact*(-1)
    504                  endif
    505                  sp_count_r=sp_count_r+1
    506                  sparse_dump_r(sp_count_r)= &
    507                       sp_fact* &
    508                       1.e12*drygrid(ix,jy)/area(ix,jy)
    509   !                sparse_dump_u(sp_count_r)=
    510   !    +            ,1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
    511               else ! concentration is zero
     523                  endif
     524                  sp_count_r=sp_count_r+1
     525                  sparse_dump_r(sp_count_r)= &
     526                       sp_fact* &
     527                       1.e12*drygrid(ix,jy)/area(ix,jy)
     528!                sparse_dump_u(sp_count_r)=
     529!    +            ,1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
     530                else ! concentration is zero
    512531                  sp_zer=.true.
    513               endif
    514             end do
    515           end do
    516          else
    517            sp_count_i=0
    518            sp_count_r=0
    519          endif
    520          write(unitoutgridppt) sp_count_i
    521          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
    522          write(unitoutgridppt) sp_count_r
    523          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
    524   !       write(unitoutgridppt) sp_count_u
    525   !       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
    526 
    527 
    528   ! Mixing ratios
    529          sp_count_i=0
    530          sp_count_r=0
    531          sp_fact=-1.
    532          sp_zer=.true.
     532                endif
     533              end do
     534            end do
     535          else
     536            sp_count_i=0
     537            sp_count_r=0
     538          endif
     539          write(unitoutgridppt) sp_count_i
     540          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
     541          write(unitoutgridppt) sp_count_r
     542          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
     543!       write(unitoutgridppt) sp_count_u
     544!       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
     545
     546
     547! Mixing ratios
     548          sp_count_i=0
     549          sp_count_r=0
     550          sp_fact=-1.
     551          sp_zer=.true.
    533552          do kz=1,numzgrid
    534553            do jy=0,numygrid-1
     
    541560                    sp_zer=.false.
    542561                    sp_fact=sp_fact*(-1.)
    543                  endif
    544                  sp_count_r=sp_count_r+1
    545                  sparse_dump_r(sp_count_r)= &
    546                       sp_fact* &
    547                       1.e12*grid(ix,jy,kz) &
    548                       /volume(ix,jy,kz)/outnum* &
    549                       weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
    550   !                sparse_dump_u(sp_count_r)=
    551   !+              ,1.e12*gridsigma(ix,jy,kz,ks,kp,nage)/volume(ix,jy,kz)/
    552   !+              outnum*weightair/weightmolar(ks)/
    553   !+              densityoutgrid(ix,jy,kz)
    554               else ! concentration is zero
     562                  endif
     563                  sp_count_r=sp_count_r+1
     564                  sparse_dump_r(sp_count_r)= &
     565                       sp_fact* &
     566                       1.e12*grid(ix,jy,kz) &
     567                       /volume(ix,jy,kz)/outnum* &
     568                       weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
     569!                sparse_dump_u(sp_count_r)=
     570!+              ,1.e12*gridsigma(ix,jy,kz,ks,kp,nage)/volume(ix,jy,kz)/
     571!+              outnum*weightair/weightmolar(ks)/
     572!+              densityoutgrid(ix,jy,kz)
     573                else ! concentration is zero
    555574                  sp_zer=.true.
    556               endif
     575                endif
    557576              end do
    558577            end do
    559578          end do
    560          write(unitoutgridppt) sp_count_i
    561          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
    562          write(unitoutgridppt) sp_count_r
    563          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
    564   !       write(unitoutgridppt) sp_count_u
    565   !       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
    566 
    567       endif ! output for ppt
    568 
    569   end do
    570   end do
     579          write(unitoutgridppt) sp_count_i
     580          write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
     581          write(unitoutgridppt) sp_count_r
     582          write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
     583!       write(unitoutgridppt) sp_count_u
     584!       write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
     585
     586        endif ! output for ppt
     587
     588      end do
     589    end do
    571590
    572591    close(unitoutgridppt)
     
    581600       drygridtotal
    582601
    583   ! Dump of receptor concentrations
    584 
    585     if (numreceptor.gt.0 .and. (iout.eq.2 .or. iout.eq.3)  ) then
    586       write(unitoutreceptppt) itime
    587       do ks=1,nspec
    588         write(unitoutreceptppt) (1.e12*creceptor(i,ks)/outnum* &
    589              weightair/weightmolar(ks)/densityoutrecept(i),i=1,numreceptor)
     602! Dump of receptor concentrations
     603
     604  if (numreceptor.gt.0 .and. (iout.eq.2 .or. iout.eq.3)  ) then
     605    write(unitoutreceptppt) itime
     606    do ks=1,nspec
     607      write(unitoutreceptppt) (1.e12*creceptor(i,ks)/outnum* &
     608           weightair/weightmolar(ks)/densityoutrecept(i),i=1,numreceptor)
     609    end do
     610  endif
     611
     612! Dump of receptor concentrations
     613
     614  if (numreceptor.gt.0) then
     615    write(unitoutrecept) itime
     616    do ks=1,nspec
     617      write(unitoutrecept) (1.e12*creceptor(i,ks)/outnum, &
     618           i=1,numreceptor)
     619    end do
     620  endif
     621
     622
     623
     624! Reinitialization of grid
     625!*************************
     626
     627  do ks=1,nspec
     628    do kp=1,maxpointspec_act
     629      do i=1,numreceptor
     630        creceptor(i,ks)=0.
    590631      end do
    591     endif
    592 
    593   ! Dump of receptor concentrations
    594 
    595     if (numreceptor.gt.0) then
    596       write(unitoutrecept) itime
    597       do ks=1,nspec
    598         write(unitoutrecept) (1.e12*creceptor(i,ks)/outnum, &
    599              i=1,numreceptor)
    600       end do
    601     endif
    602 
    603 
    604 
    605   ! Reinitialization of grid
    606   !*************************
    607 
    608   do ks=1,nspec
    609   do kp=1,maxpointspec_act
    610     do i=1,numreceptor
    611       creceptor(i,ks)=0.
    612     end do
    613     do jy=0,numygrid-1
    614       do ix=0,numxgrid-1
    615         do l=1,nclassunc
    616           do nage=1,nageclass
    617             do kz=1,numzgrid
    618               gridunc(ix,jy,kz,ks,kp,l,nage)=0.
     632      do jy=0,numygrid-1
     633        do ix=0,numxgrid-1
     634          do l=1,nclassunc
     635            do nage=1,nageclass
     636              do kz=1,numzgrid
     637                gridunc(ix,jy,kz,ks,kp,l,nage)=0.
     638              end do
    619639            end do
    620640          end do
     
    623643    end do
    624644  end do
    625   end do
    626645
    627646
  • src/concoutput_mpi.f90

    r478e9e6 rb069789  
    106106  character(len=3) :: anspec
    107107  integer :: mind
    108 ! mind        eso:added to ensure identical results between 2&3-fields versions
    109   CHARACTER(LEN=8),save :: file_stat='REPLACE'
     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
    110114
    111115! Measure execution time
    112116  if (mp_measure_time) call mpif_mtime('rootonly',0)
    113 
    114117
    115118! Determine current calendar date, needed for the file name
     
    120123  write(adate,'(i8.8)') jjjjmmdd
    121124  write(atime,'(i6.6)') ihmmss
    122   OPEN(unitdates,file=path(2)(1:length(2))//'dates', ACCESS='APPEND', STATUS=file_stat)
     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)
    123147  write(unitdates,'(a)') adate//atime
    124148  close(unitdates) 
     
    126150  ! Overwrite existing dates file on first call, later append to it
    127151  ! This fixes a bug where the dates file kept growing across multiple runs
    128   ! TODO check if the 'always APPEND'-behaviour is useful in other scenarioes
    129   ! e.g. (restart?)
    130152  IF (init) THEN
    131153    file_stat='OLD'
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG