source: trunk/src/timemanager.f90 @ 31

Last change on this file since 31 was 30, checked in by hasod, 9 years ago

ADD: Optional (compressed) netcdf output added. Activated via COMMAND
file. During compile time switches -DNETCDF_OUTPUT -cpp need to be
invoked. Compliation and linking is shown in makefile.netcdf

File size: 24.4 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 timemanager
23
24  !*****************************************************************************
25  !                                                                            *
26  ! Handles the computation of trajectories, i.e. determines which             *
27  ! trajectories have to be computed at what time.                             *
28  ! Manages dry+wet deposition routines, radioactive decay and the computation *
29  ! of concentrations.                                                         *
30  !                                                                            *
31  !     Author: A. Stohl                                                       *
32  !                                                                            *
33  !     20 May 1996                                                            *
34  !                                                                            *
35  !*****************************************************************************
36  !  Changes, Bernd C. Krueger, Feb. 2001:                                     *
37  !        Call of convmix when new windfield is read                          *
38  !------------------------------------                                        *
39  !  Changes Petra Seibert, Sept 2002                                          *
40  !     fix wet scavenging problem                                             *
41  !     Code may not be correct for decay of deposition!                       *
42  !  Changes Petra Seibert, Nov 2002                                           *
43  !     call convection BEFORE new fields are read in BWD mode                 *
44  !  Changes Caroline Forster, Feb 2005                                        *
45  !   new interface between flexpart and convection scheme                     *
46  !   Emanuel's latest subroutine convect43c.f is used                         *
47  !  Changes Stefan Henne, Harald Sodemann, 2013-2014                          *
48  !   added netcdf output code                                                 *
49  !*****************************************************************************
50  !                                                                            *
51  ! Variables:                                                                 *
52  ! dep                .true. if either wet or dry deposition is switched on   *
53  ! decay(maxspec) [1/s] decay constant for radioactive decay                  *
54  ! drydep             .true. if dry deposition is switched on                 *
55  ! ideltas [s]        modelling period                                        *
56  ! itime [s]          actual temporal position of calculation                 *
57  ! ldeltat [s]        time since computation of radioact. decay of depositions*
58  ! loutaver [s]       averaging period for concentration calculations         *
59  ! loutend [s]        end of averaging for concentration calculations         *
60  ! loutnext [s]       next time at which output fields shall be centered      *
61  ! loutsample [s]     sampling interval for averaging of concentrations       *
62  ! loutstart [s]      start of averaging for concentration calculations       *
63  ! loutstep [s]       time interval for which concentrations shall be         *
64  !                    calculated                                              *
65  ! npoint(maxpart)    index, which starting point the trajectory has          *
66  !                    starting positions of trajectories                      *
67  ! nstop              serves as indicator for fate of particles               *
68  !                    in the particle loop                                    *
69  ! nstop1             serves as indicator for wind fields (see getfields)     *
70  ! outnum             number of samples for each concentration calculation    *
71  ! outnum             number of samples for each concentration calculation    *
72  ! prob               probability of absorption at ground due to dry          *
73  !                    deposition                                              *
74  ! wetdep             .true. if wet deposition is switched on                 *
75  ! weight             weight for each concentration sample (1/2 or 1)         *
76  ! uap(maxpart),ucp(maxpart),uzp(maxpart) = random velocities due to          *
77  !                    turbulence                                              *
78  ! us(maxpart),vs(maxpart),ws(maxpart) = random velocities due to inter-      *
79  !                    polation                                                *
80  ! xtra1(maxpart), ytra1(maxpart), ztra1(maxpart) =                           *
81  !                    spatial positions of trajectories                       *
82  !                                                                            *
83  ! Constants:                                                                 *
84  ! maxpart            maximum number of trajectories                          *
85  !                                                                            *
86  !*****************************************************************************
87
88  use unc_mod
89  use point_mod
90  use xmass_mod
91  use flux_mod
92  use outg_mod
93  use oh_mod
94  use par_mod
95  use com_mod
96#ifdef NETCDF_OUTPUT
97  use netcdf_output_mod, only: concoutput_netcdf, concoutput_nest_netcdf,concoutput_surf_netcdf, concoutput_surf_nest_netcdf
98#endif
99
100  implicit none
101
102  integer :: j,ks,kp,l,n,itime,nstop,nstop1
103! integer :: ksp
104  integer :: loutnext,loutstart,loutend
105  integer :: ix,jy,ldeltat,itage,nage
106  real :: outnum,weight,prob(maxspec)
107  real :: uap(maxpart),ucp(maxpart),uzp(maxpart),decfact
108  real :: us(maxpart),vs(maxpart),ws(maxpart)
109  integer(kind=2) :: cbt(maxpart)
110  real :: drydeposit(maxspec),gridtotalunc,wetgridtotalunc
111  real :: drygridtotalunc,xold,yold,zold,xmassfract
112  !double precision xm(maxspec,maxpointspec_act),
113  !    +                 xm_depw(maxspec,maxpointspec_act),
114  !    +                 xm_depd(maxspec,maxpointspec_act)
115
116
117  !open(88,file='TEST.dat')
118
119  ! First output for time 0
120  !************************
121
122  loutnext=loutstep/2
123  outnum=0.
124  loutstart=loutnext-loutaver/2
125  loutend=loutnext+loutaver/2
126
127  !  open(127,file=path(2)(1:length(2))//'depostat.dat'
128  !    +  ,form='unformatted')
129  !write (*,*) 'writing deposition statistics depostat.dat!'
130
131  !**********************************************************************
132  ! Loop over the whole modelling period in time steps of mintime seconds
133  !**********************************************************************
134
135
136  itime=0
137  !write(*,45) itime,numpart,gridtotalunc,wetgridtotalunc,drygridtotalunc
138  write(*,46) float(itime)/3600,itime,numpart
139  if (verbosity.gt.0) then
140    write (*,*) 'timemanager> starting simulation'
141    if (verbosity.gt.1) then
142      CALL SYSTEM_CLOCK(count_clock)
143      WRITE(*,*) 'timemanager> SYSTEM CLOCK',(count_clock - count_clock0)/real(count_rate)
144    endif     
145  endif
146
147  do itime=0,ideltas,lsynctime
148
149  ! Computation of wet deposition, OH reaction and mass transfer
150  ! between two species every lsynctime seconds
151  ! maybe wet depo frequency can be relaxed later but better be on safe side
152  ! wetdepo must be called BEFORE new fields are read in but should not
153  ! be called in the very beginning before any fields are loaded, or
154  ! before particles are in the system
155  ! Code may not be correct for decay of deposition
156  ! changed by Petra Seibert 9/02
157  !********************************************************************
158
159    if (wetdep .and. itime .ne. 0 .and. numpart .gt. 0) then
160        if (verbosity.gt.0) then
161           write (*,*) 'timemanager> call wetdepo'
162        endif     
163         call wetdepo(itime,lsynctime,loutnext)
164    endif
165
166    if (ohrea .and. itime .ne. 0 .and. numpart .gt. 0) &
167         call ohreaction(itime,lsynctime,loutnext)
168
169    if (assspec .and. itime .ne. 0 .and. numpart .gt. 0) then
170       stop 'associated species not yet implemented!'
171  !     call transferspec(itime,lsynctime,loutnext)
172    endif
173
174  ! compute convection for backward runs
175  !*************************************
176
177   if ((ldirect.eq.-1).and.(lconvection.eq.1).and.(itime.lt.0)) then
178        if (verbosity.gt.0) then
179           write (*,*) 'timemanager> call convmix -- backward'
180        endif         
181      call convmix(itime)
182        if (verbosity.gt.1) then
183          !CALL SYSTEM_CLOCK(count_clock, count_rate, count_max)
184          CALL SYSTEM_CLOCK(count_clock)
185          WRITE(*,*) 'timemanager> SYSTEM CLOCK',(count_clock - count_clock0)/real(count_rate)
186        endif
187   endif
188
189  ! Get necessary wind fields if not available
190  !*******************************************
191    if (verbosity.gt.0) then
192           write (*,*) 'timemanager> call getfields'
193    endif
194    call getfields(itime,nstop1)
195        if (verbosity.gt.1) then
196          CALL SYSTEM_CLOCK(count_clock)
197          WRITE(*,*) 'timemanager> SYSTEM CLOCK',(count_clock - count_clock0)/real(count_rate)
198        endif
199    if (nstop1.gt.1) stop 'NO METEO FIELDS AVAILABLE'
200
201  ! Release particles
202  !******************
203
204    if (verbosity.gt.0) then
205           write (*,*) 'timemanager>  Release particles'
206    endif
207
208    if (mdomainfill.ge.1) then
209      if (itime.eq.0) then
210        if (verbosity.gt.0) then
211          write (*,*) 'timemanager>  call init_domainfill'
212        endif       
213        call init_domainfill
214      else
215        if (verbosity.gt.0) then
216          write (*,*) 'timemanager>  call boundcond_domainfill'
217        endif   
218        call boundcond_domainfill(itime,loutend)
219      endif
220    else
221      if (verbosity.gt.0) then
222        print*,'call releaseparticles' 
223      endif
224      call releaseparticles(itime)
225      if (verbosity.gt.1) then
226        CALL SYSTEM_CLOCK(count_clock)
227        WRITE(*,*) 'timemanager> SYSTEM CLOCK',(count_clock - count_clock0)/real(count_rate)
228      endif
229    endif
230
231
232  ! Compute convective mixing for forward runs
233  ! for backward runs it is done before next windfield is read in
234  !**************************************************************
235
236   if ((ldirect.eq.1).and.(lconvection.eq.1)) then
237     if (verbosity.gt.0) then
238       write (*,*) 'timemanager> call convmix -- forward'
239     endif   
240     call convmix(itime)
241   endif
242
243  ! If middle of averaging period of output fields is reached, accumulated
244  ! deposited mass radioactively decays
245  !***********************************************************************
246
247    if (dep.and.(itime.eq.loutnext).and.(ldirect.gt.0)) then
248      do ks=1,nspec
249      do kp=1,maxpointspec_act
250        if (decay(ks).gt.0.) then
251          do nage=1,nageclass
252            do l=1,nclassunc
253  ! Mother output grid
254              do jy=0,numygrid-1
255                do ix=0,numxgrid-1
256                  wetgridunc(ix,jy,ks,kp,l,nage)= &
257                       wetgridunc(ix,jy,ks,kp,l,nage)* &
258                       exp(-1.*outstep*decay(ks))
259                  drygridunc(ix,jy,ks,kp,l,nage)= &
260                       drygridunc(ix,jy,ks,kp,l,nage)* &
261                       exp(-1.*outstep*decay(ks))
262                end do
263              end do
264  ! Nested output grid
265              if (nested_output.eq.1) then
266                do jy=0,numygridn-1
267                  do ix=0,numxgridn-1
268                    wetgriduncn(ix,jy,ks,kp,l,nage)= &
269                         wetgriduncn(ix,jy,ks,kp,l,nage)* &
270                         exp(-1.*outstep*decay(ks))
271                    drygriduncn(ix,jy,ks,kp,l,nage)= &
272                         drygriduncn(ix,jy,ks,kp,l,nage)* &
273                         exp(-1.*outstep*decay(ks))
274                  end do
275                end do
276              endif
277            end do
278          end do
279        endif
280      end do
281      end do
282    endif
283
284  !!! CHANGE: These lines may be switched on to check the conservation
285  !!! of mass within FLEXPART
286  !   if (itime.eq.loutnext) then
287  !   do 247 ksp=1, nspec
288  !   do 247 kp=1, maxpointspec_act
289  !47         xm(ksp,kp)=0.
290
291  !   do 249 ksp=1, nspec
292  !     do 249 j=1,numpart
293  !          if (ioutputforeachrelease.eq.1) then
294  !            kp=npoint(j)
295  !          else
296  !            kp=1
297  !          endif
298  !       if (itra1(j).eq.itime) then
299  !          xm(ksp,kp)=xm(ksp,kp)+xmass1(j,ksp)
300  !         write(*,*) 'xmass: ',xmass1(j,ksp),j,ksp,nspec
301  !       endif
302  !49     continue
303  !  do 248 ksp=1,nspec
304  !  do 248 kp=1,maxpointspec_act
305  !  xm_depw(ksp,kp)=0.
306  !  xm_depd(ksp,kp)=0.
307  !     do 248 nage=1,nageclass
308  !       do 248 ix=0,numxgrid-1
309  !         do 248 jy=0,numygrid-1
310  !           do 248 l=1,nclassunc
311  !              xm_depw(ksp,kp)=xm_depw(ksp,kp)
312  !    +                  +wetgridunc(ix,jy,ksp,kp,l,nage)
313  !48                 xm_depd(ksp,kp)=xm_depd(ksp,kp)
314  !    +                  +drygridunc(ix,jy,ksp,kp,l,nage)
315  !             do 246 ksp=1,nspec
316  !46                    write(88,'(2i10,3e12.3)')
317  !    +              itime,ksp,(xm(ksp,kp),kp=1,maxpointspec_act),
318  !    +                (xm_depw(ksp,kp),kp=1,maxpointspec_act),
319  !    +                (xm_depd(ksp,kp),kp=1,maxpointspec_act)
320  !  endif
321  !!! CHANGE
322
323
324
325  ! Check whether concentrations are to be calculated
326  !**************************************************
327
328    if ((ldirect*itime.ge.ldirect*loutstart).and. &
329         (ldirect*itime.le.ldirect*loutend)) then ! add to grid
330      if (mod(itime-loutstart,loutsample).eq.0) then
331
332  ! If we are exactly at the start or end of the concentration averaging interval,
333  ! give only half the weight to this sample
334  !*****************************************************************************
335
336        if ((itime.eq.loutstart).or.(itime.eq.loutend)) then
337          weight=0.5
338        else
339          weight=1.0
340        endif
341        outnum=outnum+weight
342        call conccalc(itime,weight)
343      endif
344
345
346      if ((mquasilag.eq.1).and.(itime.eq.(loutstart+loutend)/2)) &
347           call partoutput_short(itime)    ! dump particle positions in extremely compressed format
348
349
350  ! Output and reinitialization of grid
351  ! If necessary, first sample of new grid is also taken
352  !*****************************************************
353
354      if ((itime.eq.loutend).and.(outnum.gt.0.)) then
355        if ((iout.le.3.).or.(iout.eq.5)) then
356          if (surf_only.ne.1) then
357            if (lnetcdfout.eq.1) then
358#ifdef NETCDF_OUTPUT
359              call concoutput_netcdf(itime,outnum,gridtotalunc,wetgridtotalunc,drygridtotalunc)
360#endif
361            else
362              call concoutput(itime,outnum,gridtotalunc,wetgridtotalunc,drygridtotalunc)
363            endif
364          else 
365            if (verbosity.eq.1) then
366             print*,'call concoutput_surf '
367             call system_clock(count_clock)
368             write(*,*) 'system clock',count_clock - count_clock0   
369            endif
370            if (lnetcdfout.eq.1) then
371#ifdef NETCDF_OUTPUT
372              call concoutput_surf_netcdf(itime,outnum,gridtotalunc,wetgridtotalunc,drygridtotalunc)
373#endif
374            else
375              call concoutput_surf(itime,outnum,gridtotalunc,wetgridtotalunc,drygridtotalunc)
376              if (verbosity.eq.1) then
377                print*,'called concoutput_surf '
378                call system_clock(count_clock)
379                write(*,*) 'system clock',count_clock - count_clock0   
380              endif
381            endif
382          endif
383
384          if (nested_output .eq. 1) then
385            if (lnetcdfout.eq.0) then
386              if (surf_only.ne.1) then
387                call concoutput_nest(itime,outnum)
388              else
389                call concoutput_surf_nest(itime,outnum)
390              endif
391            else
392#ifdef NETCDF_OUTPUT
393              if (surf_only.ne.1) then
394                call concoutput_nest_netcdf(itime,outnum)
395              else
396                call concoutput_surf_nest_netcdf(itime,outnum)
397              endif
398#endif
399            endif
400          endif
401          outnum=0.
402        endif
403        if ((iout.eq.4).or.(iout.eq.5)) call plumetraj(itime)
404        if (iflux.eq.1) call fluxoutput(itime)
405        !write(*,45) itime,numpart,gridtotalunc,wetgridtotalunc,drygridtotalunc
406        write(*,46) float(itime)/3600,itime,numpart
40745      format(i9,' SECONDS SIMULATED: ',i8, ' PARTICLES:    Uncertainty: ',3f7.3)
40846      format(' Simulated ',f7.1,' hours (',i9,' s), ',i8, ' particles')
409        if (ipout.ge.1) call partoutput(itime)    ! dump particle positions
410        loutnext=loutnext+loutstep
411        loutstart=loutnext-loutaver/2
412        loutend=loutnext+loutaver/2
413        if (itime.eq.loutstart) then
414          weight=0.5
415          outnum=outnum+weight
416          call conccalc(itime,weight)
417        endif
418
419
420  ! Check, whether particles are to be split:
421  ! If so, create new particles and attribute all information from the old
422  ! particles also to the new ones; old and new particles both get half the
423  ! mass of the old ones
424  !************************************************************************
425
426        if (ldirect*itime.ge.ldirect*itsplit) then
427          n=numpart
428          do j=1,numpart
429            if (ldirect*itime.ge.ldirect*itrasplit(j)) then
430              if (n.lt.maxpart) then
431                n=n+1
432                itrasplit(j)=2*(itrasplit(j)-itramem(j))+itramem(j)
433                itrasplit(n)=itrasplit(j)
434                itramem(n)=itramem(j)
435                itra1(n)=itra1(j)
436                idt(n)=idt(j)
437                npoint(n)=npoint(j)
438                nclass(n)=nclass(j)
439                xtra1(n)=xtra1(j)
440                ytra1(n)=ytra1(j)
441                ztra1(n)=ztra1(j)
442                uap(n)=uap(j)
443                ucp(n)=ucp(j)
444                uzp(n)=uzp(j)
445                us(n)=us(j)
446                vs(n)=vs(j)
447                ws(n)=ws(j)
448                cbt(n)=cbt(j)
449                do ks=1,nspec
450                  xmass1(j,ks)=xmass1(j,ks)/2.
451                  xmass1(n,ks)=xmass1(j,ks)
452                end do
453              endif
454            endif
455          end do
456          numpart=n
457        endif
458      endif
459    endif
460
461
462    if (itime.eq.ideltas) exit         ! almost finished
463
464  ! Compute interval since radioactive decay of deposited mass was computed
465  !************************************************************************
466
467    if (itime.lt.loutnext) then
468      ldeltat=itime-(loutnext-loutstep)
469    else                                  ! first half of next interval
470      ldeltat=itime-loutnext
471    endif
472
473
474  ! Loop over all particles
475  !************************
476
477    do j=1,numpart
478
479
480  ! If integration step is due, do it
481  !**********************************
482
483      if (itra1(j).eq.itime) then
484
485        if (ioutputforeachrelease.eq.1) then
486            kp=npoint(j)
487        else
488            kp=1
489        endif
490  ! Determine age class of the particle
491        itage=abs(itra1(j)-itramem(j))
492        do nage=1,nageclass
493          if (itage.lt.lage(nage)) exit
494        end do
495
496  ! Initialize newly released particle
497  !***********************************
498
499        if ((itramem(j).eq.itime).or.(itime.eq.0)) &
500             call initialize(itime,idt(j),uap(j),ucp(j),uzp(j), &
501             us(j),vs(j),ws(j),xtra1(j),ytra1(j),ztra1(j),cbt(j))
502
503  ! Memorize particle positions
504  !****************************
505
506        xold=real(xtra1(j))
507        yold=real(ytra1(j))
508        zold=ztra1(j)
509
510  ! Integrate Lagevin equation for lsynctime seconds
511  !*************************************************
512
513        call advance(itime,npoint(j),idt(j),uap(j),ucp(j),uzp(j), &
514             us(j),vs(j),ws(j),nstop,xtra1(j),ytra1(j),ztra1(j),prob, &
515             cbt(j))
516
517  ! Calculate the gross fluxes across layer interfaces
518  !***************************************************
519
520        if (iflux.eq.1) call calcfluxes(nage,j,xold,yold,zold)
521
522
523  ! Determine, when next time step is due
524  ! If trajectory is terminated, mark it
525  !**************************************
526
527        if (nstop.gt.1) then
528          if (linit_cond.ge.1) call initial_cond_calc(itime,j)
529          itra1(j)=-999999999
530        else
531          itra1(j)=itime+lsynctime
532
533
534  ! Dry deposition and radioactive decay for each species
535  ! Also check maximum (of all species) of initial mass remaining on the particle;
536  ! if it is below a threshold value, terminate particle
537  !*****************************************************************************
538
539          xmassfract=0.
540          do ks=1,nspec
541            if (decay(ks).gt.0.) then             ! radioactive decay
542              decfact=exp(-real(abs(lsynctime))*decay(ks))
543            else
544              decfact=1.
545            endif
546
547            if (drydepspec(ks)) then        ! dry deposition
548              drydeposit(ks)=xmass1(j,ks)*prob(ks)*decfact
549              xmass1(j,ks)=xmass1(j,ks)*(1.-prob(ks))*decfact
550              if (decay(ks).gt.0.) then   ! correct for decay (see wetdepo)
551                drydeposit(ks)=drydeposit(ks)* &
552                     exp(real(abs(ldeltat))*decay(ks))
553              endif
554            else                           ! no dry deposition
555              xmass1(j,ks)=xmass1(j,ks)*decfact
556            endif
557
558            if (mdomainfill.eq.0) then
559              if (xmass(npoint(j),ks).gt.0.) &
560                   xmassfract=max(xmassfract,real(npart(npoint(j)))* &
561                   xmass1(j,ks)/xmass(npoint(j),ks))
562            else
563              xmassfract=1.
564            endif
565          end do
566
567          if (xmassfract.lt.0.0001) then   ! terminate all particles carrying less mass
568            itra1(j)=-999999999
569            if (verbosity.gt.0) then
570              print*,'terminated particle ',j,' for small mass'
571            endif
572          endif
573
574  !        Sabine Eckhardt, June 2008
575  !        don't create depofield for backward runs
576          if (drydep.AND.(ldirect.eq.1)) then
577            call drydepokernel(nclass(j),drydeposit,real(xtra1(j)),real(ytra1(j)),nage,kp)
578            if (nested_output.eq.1) then
579              call drydepokernel_nest(nclass(j),drydeposit,real(xtra1(j)),real(ytra1(j)),nage,kp)
580            endif
581          endif
582
583  ! Terminate trajectories that are older than maximum allowed age
584  !***************************************************************
585
586          if (abs(itra1(j)-itramem(j)).ge.lage(nageclass)) then
587            if (linit_cond.ge.1) call initial_cond_calc(itime+lsynctime,j)
588            itra1(j)=-999999999
589            if (verbosity.gt.0) then
590              print*,'terminated particle ',j,' for age'
591            endif
592          endif
593        endif
594
595      endif
596
597    end do
598
599  end do
600
601
602  ! Complete the calculation of initial conditions for particles not yet terminated
603  !*****************************************************************************
604
605  do j=1,numpart
606    if (linit_cond.ge.1) call initial_cond_calc(itime,j)
607  end do
608
609  if (ipout.eq.2) call partoutput(itime)     ! dump particle positions
610
611  if (linit_cond.ge.1) call initial_cond_output(itime)   ! dump initial cond. field
612
613  close(104)
614
615  ! De-allocate memory and end
616  !***************************
617
618  if (iflux.eq.1) then
619    deallocate(flux)
620  endif
621  if (ohrea.eqv..TRUE.) then
622    deallocate(OH_field,OH_field_height)
623  endif
624  if (ldirect.gt.0) then
625    deallocate(drygridunc,wetgridunc)
626  endif
627  deallocate(gridunc)
628  deallocate(xpoint1,xpoint2,ypoint1,ypoint2,zpoint1,zpoint2,xmass)
629  deallocate(ireleasestart,ireleaseend,npart,kindz)
630  deallocate(xmasssave)
631  if (nested_output.eq.1) then
632    deallocate(orooutn, arean, volumen)
633    if (ldirect.gt.0) then
634      deallocate(griduncn,drygriduncn,wetgriduncn)
635    endif
636  endif
637  deallocate(outheight,outheighthalf)
638  deallocate(oroout,area,volume)
639
640end subroutine timemanager
641
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG