source: trunk/src/partdep.f90 @ 28

Last change on this file since 28 was 4, checked in by mlanger, 11 years ago
  • Property svn:executable set to *
File size: 6.5 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 partdep(nc,density,fract,schmi,vset,ra,ustar,nyl,vdep)
23  !                   i     i      i     i    i   i    i    i  i/o
24  !*****************************************************************************
25  !                                                                            *
26  !      Calculation of the dry deposition velocities of particles.            *
27  !      This routine is based on Stokes' law for considering settling and     *
28  !      assumes constant dynamic viscosity of the air.                        *
29  !                                                                            *
30  !     AUTHOR: Andreas Stohl, 12 November 1993                                *
31  !                            Update: 20 December 1996                        *
32  !                                                                            *
33  !     Literature:                                                            *
34  !     [1]  Hicks/Baldocchi/Meyers/Hosker/Matt (1987), A Preliminary          *
35  !             Multiple Resistance Routine for Deriving Dry Deposition        *
36  !             Velocities from Measured Quantities.                           *
37  !             Water, Air and Soil Pollution 36 (1987), pp.311-330.           *
38  !     [2]  Slinn (1982), Predictions for Particle Deposition to              *
39  !             Vegetative Canopies. Atm.Env.16-7 (1982), pp.1785-1794.        *
40  !     [3]  Slinn/Slinn (1980),  Predictions for Particle Deposition on       *
41  !             Natural Waters. Atm.Env.14 (1980), pp.1013-1016.               *
42  !     [4]  Scire/Yamartino/Carmichael/Chang (1989),                          *
43  !             CALGRID: A Mesoscale Photochemical Grid Model.                 *
44  !             Vol II: User's Guide. (Report No.A049-1, June, 1989)           *
45  !     [5]  Langer M. (1992): Ein einfaches Modell zur Abschaetzung der       *
46  !             Depositionsgeschwindigkeit von Teilchen und Gasen.             *
47  !             Internal report.                                               *
48  !                                                                            *
49  !*****************************************************************************
50  !                                                                            *
51  ! Variables:                                                                 *
52  ! alpha                help variable                                         *
53  ! fract(nc,ni)         mass fraction of each diameter interval               *
54  ! lpdep(nc)            1 for particle deposition, 0 else                     *
55  ! nc                   actual number of chemical components                  *
56  ! ni                   number of diameter intervals, for which vdepj is calc.*
57  ! rdp [s/m]            deposition layer resistance                           *
58  ! ra [s/m]             aerodynamical resistance                              *
59  ! schmi(nc,ni)         Schmidt number**2/3 of each diameter interval         *
60  ! stokes               Stokes number                                         *
61  ! ustar [m/s]          friction velocity                                     *
62  ! vdep(nc) [m/s]       deposition velocities of all components               *
63  ! vdepj [m/s]          help, deposition velocity of 1 interval               *
64  ! vset(nc,ni)          gravitational settling velocity of each interval      *
65  !                                                                            *
66  ! Constants:                                                                 *
67  ! nc                   number of chemical species                            *
68  ! ni                   number of diameter intervals, for which deposition    *
69  !                      is calculated                                         *
70  !                                                                            *
71  !*****************************************************************************
72
73  use par_mod
74
75  implicit none
76
77  real :: density(maxspec),schmi(maxspec,ni),fract(maxspec,ni)
78  real :: vset(maxspec,ni)
79  real :: vdep(maxspec),stokes,vdepj,rdp,ustar,alpha,ra,nyl
80  real,parameter :: eps=1.e-5
81  integer :: ic,j,nc
82
83
84  do ic=1,nc                  ! loop over all species
85    if (density(ic).gt.0.) then
86      do j=1,ni              ! loop over all diameter intervals
87        if (ustar.gt.eps) then
88
89  ! Stokes number for each diameter interval
90  !*****************************************
91
92          stokes=vset(ic,j)/ga*ustar*ustar/nyl
93          alpha=-3./stokes
94
95  ! Deposition layer resistance
96  !****************************
97
98          if (alpha.le.log10(eps)) then
99            rdp=1./(schmi(ic,j)*ustar)
100          else
101                rdp=1./((schmi(ic,j)+10.**alpha)*ustar)
102          endif
103          vdepj=vset(ic,j)+1./(ra+rdp+ra*rdp*vset(ic,j))
104        else
105          vdepj=vset(ic,j)
106        endif
107
108  ! deposition velocities of each interval are weighted with mass fraction
109  !***********************************************************************
110
111        vdep(ic)=vdep(ic)+vdepj*fract(ic,j)
112      end do
113    endif
114  end do
115
116end subroutine partdep
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG