source: trunk/src/getrb.f90 @ 27

Last change on this file since 27 was 4, checked in by mlanger, 11 years ago
File size: 3.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 getrb(nc,ustar,nyl,diffh2o,reldiff,rb)
23  !                 i    i    i     i       i    o
24  !*****************************************************************************
25  !                                                                            *
26  !  Calculation of the quasilaminar sublayer resistance to dry deposition.    *
27  !                                                                            *
28  !      AUTHOR: Andreas Stohl, 20 May 1995                                    *
29  !                                                                            *
30  !*****************************************************************************
31  !                                                                            *
32  ! Variables:                                                                 *
33  ! rb(ncmax)       sublayer resistance                                        *
34  ! schmidt         Schmidt number                                             *
35  ! ustar [m/s]     friction velocity                                          *
36  ! diffh20 [m2/s]  diffusivity of water vapor in air                          *
37  ! reldiff         diffusivity relative to H2O                                *
38  !                                                                            *
39  ! Constants:                                                                 *
40  ! karman          von Karman constant                                        *
41  ! pr              Prandtl number                                             *
42  !                                                                            *
43  !*****************************************************************************
44
45  use par_mod
46
47  implicit none
48
49  real :: ustar,diffh2o,rb(maxspec),schmidt,nyl
50  real :: reldiff(maxspec)
51  integer :: ic,nc
52  real,parameter :: pr=0.72
53
54  do ic=1,nc
55    if (reldiff(ic).gt.0.) then
56      schmidt=nyl/diffh2o*reldiff(ic)
57      rb(ic)=2.0*(schmidt/pr)**0.67/(karman*ustar)
58    endif
59  end do
60
61end subroutine getrb
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG