source: trunk/src/readreleases.f90 @ 28

Last change on this file since 28 was 27, checked in by hasod, 10 years ago
  • Implemented optional namelist input for COMMAND, RELEASES, SPECIES, AGECLASSES,OUTGRID,OUTGRID_NEST,RECEPTORS
  • Implemented com_mod switch nmlout to write input files as namelist to the output directory (.true. by default)
  • Proposed updated startup and runtime output (may change back to previous info if desired)
  • Property svn:executable set to *
File size: 24.2 KB
RevLine 
[4]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 readreleases
23
24  !*****************************************************************************
25  !                                                                            *
26  !     This routine reads the release point specifications for the current    *
27  !     model run. Several release points can be used at the same time.        *
28  !                                                                            *
29  !     Author: A. Stohl                                                       *
30  !                                                                            *
31  !     18 May 1996                                                            *
32  !                                                                            *
33  !     Update: 29 January 2001                                                *
34  !     Release altitude can be either in magl or masl                         *
[27]35  !     HSO, 12 August 2013
36  !     Added optional namelist input
[4]37  !                                                                            *
38  !*****************************************************************************
39  !                                                                            *
40  ! Variables:                                                                 *
41  ! decay               decay constant of species                              *
42  ! dquer [um]          mean particle diameters                                *
43  ! dsigma              e.g. dsigma=10 or dsigma=0.1 means that 68% of the mass*
44  !                     are between 0.1*dquer and 10*dquer                     *
45  ! ireleasestart, ireleaseend [s] starting time and ending time of each       *
46  !                     release                                                *
47  ! kindz               1: zpoint is in m agl, 2: zpoint is in m asl, 3: zpoint*
48  !                     is in hPa                                              *
49  ! npart               number of particles to be released                     *
50  ! nspec               number of species to be released                       *
51  ! density [kg/m3]     density of the particles                               *
52  ! rm [s/m]            Mesophyll resistance                                   *
53  ! species             name of species                                        *
54  ! xmass               total mass of each species                             *
55  ! xpoint1,ypoint1     geograf. coordinates of lower left corner of release   *
56  !                     area                                                   *
57  ! xpoint2,ypoint2     geograf. coordinates of upper right corner of release  *
58  !                     area                                                   *
[24]59  ! weta, wetb          parameters to determine the below-cloud scavenging     *
60  ! weta_in, wetb_in    parameters to determine the in-cloud scavenging        *
61  ! wetc_in, wetd_in    parameters to determine the in-cloud scavenging        *
[4]62  ! zpoint1,zpoint2     height range, over which release takes place           *
[20]63  ! num_min_discrete    if less, release cannot be randomized and happens at   *
64  !                     time mid-point of release interval                     *
[4]65  !                                                                            *
66  !*****************************************************************************
67
68  use point_mod
69  use xmass_mod
70  use par_mod
71  use com_mod
72
73  implicit none
74
[27]75  integer :: numpartmax,i,j,id1,it1,id2,it2,idum,stat
[20]76  integer,parameter :: num_min_discrete=100
[4]77  real :: vsh(ni),fracth(ni),schmih(ni),releaserate,xdum,cun
[20]78  real(kind=dp) :: jul1,jul2,julm,juldate
[4]79  character(len=50) :: line
80  logical :: old
81
[27]82  ! help variables for namelist reading
83  integer :: numpoints, parts, readerror
84  integer*2 :: zkind
85  integer :: idate1, itime1, idate2, itime2
86  real :: lon1,lon2,lat1,lat2,z1,z2
87  character*40 :: comment
88  integer,parameter :: unitreleasesout=2
89  real,allocatable, dimension (:) :: mass
90  integer,allocatable, dimension (:) :: specnum_rel,specnum_rel2
91
92  ! declare namelists
93  namelist /releases_ctrl/ &
94    nspec, &
95    specnum_rel
96
97  namelist /release/ &
98    idate1, itime1, &
99    idate2, itime2, &
100    lon1, lon2, &
101    lat1, lat2, &
102    z1, z2, &
103    zkind, &
104    mass, &
105    parts, &
106    comment
107
108  numpoint=0
109
110  ! allocate with maxspec for first input loop
111  allocate(mass(maxspec),stat=stat)
112  if (stat.ne.0) write(*,*)'ERROR: could not allocate mass'
113  allocate(specnum_rel(maxspec),stat=stat)
114  if (stat.ne.0) write(*,*)'ERROR: could not allocate specnum_rel'
115
116  ! presetting namelist releases_ctrl
117  nspec = -1  ! use negative value to determine failed namelist input
118  specnum_rel = 0
119
[4]120  !sec, read release to find how many releasepoints should be allocated
[27]121  open(unitreleases,file=path(1)(1:length(1))//'RELEASES',status='old',form='formatted',err=999)
[4]122
[27]123  ! check if namelist input provided
124  read(unitreleases,releases_ctrl,iostat=readerror)
[4]125
[27]126  ! prepare namelist output if requested
127  if (nmlout.eqv..true.) then
128    open(unitreleasesout,file=path(2)(1:length(2))//'RELEASES.namelist',access='append',status='new',err=1000)
[4]129  endif
130
[27]131  if ((readerror.ne.0).or.(nspec.lt.0)) then
[4]132
[27]133    ! no namelist format, close file and allow reopening in old format
134    close(unitreleases)
135    open(unitreleases,file=path(1)(1:length(1))//'RELEASES',status='old',err=999)
[4]136
[27]137    readerror=1 ! indicates old format
[4]138
[27]139    ! Check the format of the RELEASES file (either in free format,
140    ! or using a formatted mask)
141    ! Use of formatted mask is assumed if line 10 contains the word 'DIRECTION'
142    !**************************************************************************
[4]143
[27]144    call skplin(12,unitreleases)
145    read (unitreleases,901) line
146901 format (a)
147    if (index(line,'Total') .eq. 0) then
148      old = .false.
149    else
150      old = .true.
151    endif
152    rewind(unitreleases)
153
154
155    ! Skip first 11 lines (file header)
156    !**********************************
157
158    call skplin(11,unitreleases)
159
160    read(unitreleases,*,err=998) nspec
[4]161    if (old) call skplin(2,unitreleases)
[27]162    do i=1,nspec
163      read(unitreleases,*,err=998) specnum_rel(i)
164      if (old) call skplin(2,unitreleases)
165    end do
[4]166
[27]167    numpoint=0
168100 numpoint=numpoint+1
169    read(unitreleases,*,end=25)
170    read(unitreleases,*,err=998,end=25) idum,idum
171    if (old) call skplin(2,unitreleases)
172    read(unitreleases,*,err=998) idum,idum
173    if (old) call skplin(2,unitreleases)
[4]174    read(unitreleases,*,err=998) xdum
175    if (old) call skplin(2,unitreleases)
[27]176    read(unitreleases,*,err=998) xdum
177    if (old) call skplin(2,unitreleases)
178    read(unitreleases,*,err=998) xdum
179    if (old) call skplin(2,unitreleases)
180    read(unitreleases,*,err=998) xdum
181    if (old) call skplin(2,unitreleases)
182    read(unitreleases,*,err=998) idum
183    if (old) call skplin(2,unitreleases)
184    read(unitreleases,*,err=998) xdum
185    if (old) call skplin(2,unitreleases)
186    read(unitreleases,*,err=998) xdum
187    if (old) call skplin(2,unitreleases)
188    read(unitreleases,*,err=998) idum
189    if (old) call skplin(2,unitreleases)
190    do i=1,nspec
191      read(unitreleases,*,err=998) xdum
192      if (old) call skplin(2,unitreleases)
193    end do
194    !save compoint only for the first 1000 release points
195    read(unitreleases,'(a40)',err=998) compoint(1)(1:40)
196    if (old) call skplin(1,unitreleases)
[4]197
[27]198    goto 100
[4]199
[27]20025  numpoint=numpoint-1
[4]201
[27]202  else
[4]203
[27]204    readerror=0
205    do while (readerror.eq.0)
206      idate1=-1
207      read(unitreleases,release,iostat=readerror)
208      if ((idate1.lt.0).or.(readerror.ne.0)) then
209        readerror=1
210      else
211        numpoint=numpoint+1
212      endif
213    end do
214    readerror=0
215  endif ! if namelist input
[4]216
[27]217  rewind(unitreleases)
[4]218
[27]219  ! allocate arrays of matching size for number of species (namelist output)
220  deallocate(mass)
221  allocate(mass(nspec),stat=stat)
222  if (stat.ne.0) write(*,*)'ERROR: could not allocate mass'
223  allocate(specnum_rel2(nspec),stat=stat)
224  if (stat.ne.0) write(*,*)'ERROR: could not allocate specnum_rel2'
225  specnum_rel2=specnum_rel(1:nspec)
226  deallocate(specnum_rel)
227  allocate(specnum_rel(nspec),stat=stat)
228  if (stat.ne.0) write(*,*)'ERROR: could not allocate specnum_rel'
229  specnum_rel=specnum_rel2
230  deallocate(specnum_rel2)
231
232  !allocate memory for numpoint releaspoints
233  allocate(ireleasestart(numpoint),stat=stat)
234  if (stat.ne.0) write(*,*)'ERROR: could not allocate ireleasestart'
235  allocate(ireleaseend(numpoint),stat=stat)
236  if (stat.ne.0) write(*,*)'ERROR: could not allocate ireleaseend'
237  allocate(xpoint1(numpoint),stat=stat)
238  if (stat.ne.0) write(*,*)'ERROR: could not allocate xpoint1'
239  allocate(xpoint2(numpoint),stat=stat)
240  if (stat.ne.0) write(*,*)'ERROR: could not allocate xpoint2'
241  allocate(ypoint1(numpoint),stat=stat)
242  if (stat.ne.0) write(*,*)'ERROR: could not allocate ypoint1'
243  allocate(ypoint2(numpoint),stat=stat)
244  if (stat.ne.0) write(*,*)'ERROR: could not allocate ypoint2'
245  allocate(zpoint1(numpoint),stat=stat)
246  if (stat.ne.0) write(*,*)'ERROR: could not allocate zpoint1'
247  allocate(zpoint2(numpoint),stat=stat)
248  if (stat.ne.0) write(*,*)'ERROR: could not allocate zpoint2'
249  allocate(kindz(numpoint),stat=stat)
250  if (stat.ne.0) write(*,*)'ERROR: could not allocate kindz'
251  allocate(xmass(numpoint,maxspec),stat=stat)
252  if (stat.ne.0) write(*,*)'ERROR: could not allocate xmass'
253  allocate(rho_rel(numpoint),stat=stat)
254  if (stat.ne.0) write(*,*)'ERROR: could not allocate rho_rel'
255  allocate(npart(numpoint),stat=stat)
256  if (stat.ne.0) write(*,*)'ERROR: could not allocate npart'
257  allocate(xmasssave(numpoint),stat=stat)
258  if (stat.ne.0) write(*,*)'ERROR: could not allocate xmasssave'
259
260  write (*,*) 'Releasepoints : ', numpoint
261
262  do i=1,numpoint
263    xmasssave(i)=0.
264  end do
265
[4]266  !now save the information
267  DEP=.false.
268  DRYDEP=.false.
269  WETDEP=.false.
270  OHREA=.false.
271  do i=1,maxspec
272    DRYDEPSPEC(i)=.false.
273  end do
274
[27]275  if (readerror.ne.0) then
276    ! Skip first 11 lines (file header)
277    !**********************************
[4]278
[27]279    call skplin(11,unitreleases)
[4]280
[27]281    ! Assign species-specific parameters needed for physical processes
282    !*************************************************************************
[4]283
[27]284    read(unitreleases,*,err=998) nspec
285    if (nspec.gt.maxspec) goto 994
286    if (old) call skplin(2,unitreleases)
287  endif
[4]288
[27]289  ! namelist output
290  if (nmlout.eqv..true.) then
291    write(unitreleasesout,nml=releases_ctrl)
292  endif
[4]293
294  do i=1,nspec
[27]295    if (readerror.ne.0) then
296      read(unitreleases,*,err=998) specnum_rel(i)
297      if (old) call skplin(2,unitreleases)
298      call readspecies(specnum_rel(i),i)
299    else
300      call readspecies(specnum_rel(i),i)
301    endif
[4]302
[27]303    ! For backward runs, only 1 species is allowed
304    !*********************************************
[4]305
[27]306    !if ((ldirect.lt.0).and.(nspec.gt.1)) then
307    !write(*,*) '#####################################################'
308    !write(*,*) '#### FLEXPART MODEL SUBROUTINE READRELEASES:     ####'
309    !write(*,*) '#### FOR BACKWARD RUNS, ONLY 1 SPECIES IS ALLOWED####'
310    !write(*,*) '#####################################################'
311    !  stop
312    !endif
[4]313
[27]314    ! Molecular weight
315    !*****************
[4]316
[27]317    if (((iout.eq.2).or.(iout.eq.3)).and.(weightmolar(i).lt.0.)) then
[4]318      write(*,*) 'For mixing ratio output, valid molar weight'
319      write(*,*) 'must be specified for all simulated species.'
320      write(*,*) 'Check table SPECIES or choose concentration'
321      write(*,*) 'output instead if molar weight is not known.'
322      stop
323    endif
324
[27]325    ! Radioactive decay
326    !******************
[4]327
328    decay(i)=0.693147/decay(i) !conversion half life to decay constant
329
330
[27]331    ! Dry deposition of gases
332    !************************
[4]333
[27]334    if (reldiff(i).gt.0.) rm(i)=1./(henry(i)/3000.+100.*f0(i))    ! mesophyll resistance
[4]335
[27]336    ! Dry deposition of particles
337    !****************************
[4]338
339    vsetaver(i)=0.
340    cunningham(i)=0.
341    dquer(i)=dquer(i)*1000000.         ! Conversion m to um
[27]342    if (density(i).gt.0.) then         ! Additional parameters
343      call part0(dquer(i),dsigma(i),density(i),fracth,schmih,cun,vsh)
[4]344      do j=1,ni
345        fract(i,j)=fracth(j)
346        schmi(i,j)=schmih(j)
347        vset(i,j)=vsh(j)
348        cunningham(i)=cunningham(i)+cun*fract(i,j)
349        vsetaver(i)=vsetaver(i)-vset(i,j)*fract(i,j)
350      end do
[24]351      write(*,*) 'Average settling velocity: ',i,vsetaver(i)
[4]352    endif
353
[27]354    ! Dry deposition for constant deposition velocity
355    !************************************************
[4]356
357    dryvel(i)=dryvel(i)*0.01         ! conversion to m/s
358
[27]359    ! Check if wet deposition or OH reaction shall be calculated
360    !***********************************************************
[4]361    if (weta(i).gt.0.)  then
362      WETDEP=.true.
[27]363      write (*,*) 'Below-cloud scavenging: ON'
[24]364      write (*,*) 'Below-cloud scavenging coefficients: ',weta(i),i
365    else
[27]366      write (*,*) 'Below-cloud scavenging: OFF'
[4]367    endif
[24]368   
[27]369    ! NIK 31.01.2013 + 10.12.2013
[24]370    if (weta_in(i).gt.0.)  then
371      WETDEP=.true.
[27]372      write (*,*) 'In-cloud scavenging: ON'
[24]373      write (*,*) 'In-cloud scavenging coefficients: ',weta_in(i),wetb_in(i), wetc_in(i), wetd_in(i),i
374    else
[27]375      write (*,*) 'In-cloud scavenging: OFF'
[24]376    endif
377
[4]378    if (ohreact(i).gt.0) then
379      OHREA=.true.
[27]380      write (*,*) 'OHreaction: ON (',ohreact(i),i,')'
[4]381    endif
382
[27]383    if ((reldiff(i).gt.0.).or.(density(i).gt.0.).or.(dryvel(i).gt.0.)) then
[4]384      DRYDEP=.true.
385      DRYDEPSPEC(i)=.true.
386    endif
387
388  end do
389
[27]390  if (WETDEP.or.DRYDEP) DEP=.true.
[4]391
392  ! Read specifications for each release point
393  !*******************************************
[27]394  numpoints=numpoint
[4]395  numpoint=0
396  numpartmax=0
397  releaserate=0.
[27]398101   numpoint=numpoint+1
399
400  if (readerror.lt.1) then ! reading namelist format
401
402    if (numpoint.gt.numpoints) goto 250
403    zkind = 1
404    mass = 0
405    parts = 0
406    comment = ' '
407    read(unitreleases,release,iostat=readerror)
408    id1=idate1
409    it1=itime1
410    id2=idate2
411    it2=itime2
412    xpoint1(numpoint)=lon1
413    xpoint2(numpoint)=lon2
414    ypoint1(numpoint)=lat1
415    ypoint2(numpoint)=lat2
416    zpoint1(numpoint)=z1
417    zpoint2(numpoint)=z2
418    kindz(numpoint)=zkind
419    do i=1,nspec
420      xmass(numpoint,i)=mass(i)
421    end do
422    npart(numpoint)=parts
423    compoint(min(1001,numpoint))=comment
424
425    ! namelist output
426    if (nmlout.eqv..true.) then
427      write(unitreleasesout,nml=release)
428    endif
429
[4]430  else
431
[27]432    read(unitreleases,*,end=250)
433    read(unitreleases,*,err=998,end=250) id1,it1
434    if (old) call skplin(2,unitreleases)
435    read(unitreleases,*,err=998) id2,it2
436    if (old) call skplin(2,unitreleases)
437    read(unitreleases,*,err=998) xpoint1(numpoint)
438    if (old) call skplin(2,unitreleases)
439    read(unitreleases,*,err=998) ypoint1(numpoint)
440    if (old) call skplin(2,unitreleases)
441    read(unitreleases,*,err=998) xpoint2(numpoint)
442    if (old) call skplin(2,unitreleases)
443    read(unitreleases,*,err=998) ypoint2(numpoint)
444    if (old) call skplin(2,unitreleases)
445    read(unitreleases,*,err=998) kindz(numpoint)
446    if (old) call skplin(2,unitreleases)
447    read(unitreleases,*,err=998) zpoint1(numpoint)
448    if (old) call skplin(2,unitreleases)
449    read(unitreleases,*,err=998) zpoint2(numpoint)
450    if (old) call skplin(2,unitreleases)
451    read(unitreleases,*,err=998) npart(numpoint)
452    if (old) call skplin(2,unitreleases)
453    do i=1,nspec
454      read(unitreleases,*,err=998) xmass(numpoint,i)
455      if (old) call skplin(2,unitreleases)
456      mass(i)=xmass(numpoint,i)
457    end do
458    !save compoint only for the first 1000 release points
459    if (numpoint.le.1000) then
460      read(unitreleases,'(a40)',err=998) compoint(numpoint)(1:40)
461      comment=compoint(numpoint)(1:40)
462    else
463      read(unitreleases,'(a40)',err=998) compoint(1001)(1:40)
464      comment=compoint(1001)(1:40)
465    endif
466    if (old) call skplin(1,unitreleases)
[4]467
[27]468    ! namelist output
469    if (nmlout.eqv..true.) then
470      idate1=id1
471      itime1=it1
472      idate2=id2
473      itime2=it2
474      lon1=xpoint1(numpoint)
475      lon2=xpoint2(numpoint)
476      lat1=ypoint1(numpoint)
477      lat2=ypoint2(numpoint)
478      z1=zpoint1(numpoint)
479      z2=zpoint2(numpoint)
480      zkind=kindz(numpoint)
481      parts=npart(numpoint)
482      write(unitreleasesout,nml=release)
483    endif
484
485    if (numpoint.le.1000) then
486      if((xpoint1(numpoint).eq.0.).and.(ypoint1(numpoint).eq.0.).and. &
487           (xpoint2(numpoint).eq.0.).and.(ypoint2(numpoint).eq.0.).and. &
488           (compoint(numpoint)(1:8).eq.'        ')) goto 250
489    else
490      if((xpoint1(numpoint).eq.0.).and.(ypoint1(numpoint).eq.0.).and. &
491           (xpoint2(numpoint).eq.0.).and.(ypoint2(numpoint).eq.0.)) goto 250
492    endif
493
494  endif ! if namelist format
495
[4]496  ! If a release point contains no particles, stop and issue error message
497  !***********************************************************************
498
499  if (npart(numpoint).eq.0) then
500    write(*,*) 'FLEXPART MODEL ERROR'
501    write(*,*) 'RELEASES file is corrupt.'
502    write(*,*) 'At least for one release point, there are zero'
503    write(*,*) 'particles released. Make changes to RELEASES.'
504    stop
505  endif
506
507  ! Check whether x coordinates of release point are within model domain
508  !*********************************************************************
509
510   if (xpoint1(numpoint).lt.xlon0) &
511        xpoint1(numpoint)=xpoint1(numpoint)+360.
512   if (xpoint1(numpoint).gt.xlon0+(nxmin1)*dx) &
513        xpoint1(numpoint)=xpoint1(numpoint)-360.
514   if (xpoint2(numpoint).lt.xlon0) &
515        xpoint2(numpoint)=xpoint2(numpoint)+360.
516   if (xpoint2(numpoint).gt.xlon0+(nxmin1)*dx) &
517        xpoint2(numpoint)=xpoint2(numpoint)-360.
518
519  ! Determine relative beginning and ending times of particle release
520  !******************************************************************
521
522  jul1=juldate(id1,it1)
523  jul2=juldate(id2,it2)
[20]524  julm=(jul1+jul2)/2.
[4]525  if (jul1.gt.jul2) then
526    write(*,*) 'FLEXPART MODEL ERROR'
527    write(*,*) 'Release stops before it begins.'
528    write(*,*) 'Make changes to file RELEASES.'
529    stop
530  endif
531  if (mdomainfill.eq.0) then   ! no domain filling
532    if (ldirect.eq.1) then
533      if ((jul1.lt.bdate).or.(jul2.gt.edate)) then
534        write(*,*) 'FLEXPART MODEL ERROR'
535        write(*,*) 'Release starts before simulation begins or ends'
536        write(*,*) 'after simulation stops.'
537        write(*,*) 'Make files COMMAND and RELEASES consistent.'
538        stop
539      endif
[20]540      if (npart(numpoint).gt.num_min_discrete) then
541        ireleasestart(numpoint)=int((jul1-bdate)*86400.)
542        ireleaseend(numpoint)=int((jul2-bdate)*86400.)
543      else
544        ireleasestart(numpoint)=int((julm-bdate)*86400.)
545        ireleaseend(numpoint)=int((julm-bdate)*86400.)
546      endif
[4]547    else if (ldirect.eq.-1) then
548      if ((jul1.lt.edate).or.(jul2.gt.bdate)) then
549        write(*,*) 'FLEXPART MODEL ERROR'
550        write(*,*) 'Release starts before simulation begins or ends'
551        write(*,*) 'after simulation stops.'
552        write(*,*) 'Make files COMMAND and RELEASES consistent.'
553        stop
554      endif
[20]555      if (npart(numpoint).gt.num_min_discrete) then
556        ireleasestart(numpoint)=int((jul1-bdate)*86400.)
557        ireleaseend(numpoint)=int((jul2-bdate)*86400.)
558      else
559        ireleasestart(numpoint)=int((julm-bdate)*86400.)
560        ireleaseend(numpoint)=int((julm-bdate)*86400.)
561      endif
[4]562    endif
563  endif
564
565  ! Determine the release rate (particles per second) and total number
566  ! of particles released during the simulation
567  !*******************************************************************
568
569  if (ireleasestart(numpoint).ne.ireleaseend(numpoint)) then
570    releaserate=releaserate+real(npart(numpoint))/ &
571         real(ireleaseend(numpoint)-ireleasestart(numpoint))
572  else
573    releaserate=99999999
574  endif
575  numpartmax=numpartmax+npart(numpoint)
[27]576  goto 101
[4]577
578250   close(unitreleases)
579
[27]580  if (nmlout.eqv..true.) then
581    close(unitreleasesout)
582  endif
583
584  write (*,*) 'Particles allocated (maxpart)  : ',maxpart
585  write (*,*) 'Particles released (numpartmax): ',numpartmax
[4]586  numpoint=numpoint-1
587
588  if (ioutputforeachrelease.eq.1) then
589    maxpointspec_act=numpoint
590  else
591    maxpointspec_act=1
592  endif
593
[27]594  ! Check, whether the total number of particles may exceed totally allowed
595  ! number of particles at some time during the simulation
596  !************************************************************************
597
[4]598  if (releaserate.gt. &
599       0.99*real(maxpart)/real(lage(nageclass))) then
600    if (numpartmax.gt.maxpart) then
601  write(*,*) '#####################################################'
602  write(*,*) '#### FLEXPART MODEL SUBROUTINE READRELEASES:     ####'
603  write(*,*) '####                                             ####'
604  write(*,*) '####WARNING - TOTAL NUMBER OF PARTICLES SPECIFIED####'
605  write(*,*) '#### IN FILE "RELEASES" MAY AT SOME POINT DURING ####'
606  write(*,*) '#### THE SIMULATION EXCEED THE MAXIMUM ALLOWED   ####'
607  write(*,*) '#### NUMBER (MAXPART).IF RELEASES DO NOT OVERLAP,####'
608  write(*,*) '#### FLEXPART CAN POSSIBLY COMPLETE SUCCESSFULLY.####'
609  write(*,*) '#### HOWEVER, FLEXPART MAY HAVE TO STOP          ####'
610  write(*,*) '#### AT SOME TIME DURING THE SIMULATION. PLEASE  ####'
611  write(*,*) '#### MAKE SURE THAT YOUR SETTINGS ARE CORRECT.   ####'
612  write(*,*) '#####################################################'
613      write(*,*) 'Maximum release rate may be: ',releaserate, &
614           ' particles per second'
615      write(*,*) 'Maximum allowed release rate is: ', &
616           real(maxpart)/real(lage(nageclass)),' particles per second'
617      write(*,*) &
618           'Total number of particles released during the simulation is: ', &
619           numpartmax
620      write(*,*) 'Maximum allowed number of particles is: ',maxpart
621    endif
622  endif
623
624  return
625
626994   write(*,*) '#####################################################'
627  write(*,*) '#### FLEXPART MODEL SUBROUTINE READRELEASES:     ####'
628  write(*,*) '####                                             ####'
629  write(*,*) '#### ERROR - MAXIMUM NUMBER OF EMITTED SPECIES IS####'
630  write(*,*) '#### TOO LARGE. PLEASE REDUCE NUMBER OF SPECIES. ####'
631  write(*,*) '#####################################################'
632  stop
633
634998   write(*,*) '#####################################################'
635  write(*,*) '#### FLEXPART MODEL SUBROUTINE READRELEASES:     ####'
636  write(*,*) '####                                             ####'
637  write(*,*) '#### FATAL ERROR - FILE "RELEASES" IS            ####'
638  write(*,*) '#### CORRUPT. PLEASE CHECK YOUR INPUTS FOR       ####'
639  write(*,*) '#### MISTAKES OR GET A NEW "RELEASES"-           ####'
640  write(*,*) '#### FILE ...                                    ####'
641  write(*,*) '#####################################################'
642  stop
643
644
645999   write(*,*) '#####################################################'
646  write(*,*) '   FLEXPART MODEL SUBROUTINE READRELEASES: '
647  write(*,*)
648  write(*,*) 'FATAL ERROR - FILE CONTAINING PARTICLE RELEASE POINTS'
649  write(*,*) 'POINTS IS NOT AVAILABLE OR YOU ARE NOT'
650  write(*,*) 'PERMITTED FOR ANY ACCESS'
651  write(*,*) '#####################################################'
652  stop
653
[27]6541000 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "RELEASES"    #### '
655  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
656  write(*,'(a)') path(2)(1:length(2))
657  stop
658
[4]659end subroutine readreleases
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG