source: flexpart.git/src/getrc.f90 @ 38b7917

10.4.1_peseiFPv9.3.1FPv9.3.1b_testingFPv9.3.2GFS_025NetCDFbugfixes+enhancementsdepositiondevflexpart-noresmfp9.3.1-20161214-nc4grib2nc4_repairinputlistlaptoprelease-10release-10.4.1scaling-bugsvn-petrasvn-trunkunivie
Last change on this file since 38b7917 was e200b7a, checked in by Matthias Langer <matthias.langer@…>, 11 years ago

git-svn-id: http://flexpart.flexpart.eu:8088/svn/FlexPart90/trunk@4 ef8cc7e1-21b7-489e-abab-c1baa636049d

  • Property mode set to 100644
File size: 5.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 getrc(nc,i,j,t,gr,rh,rr,rc)
23  !                 i  i i i i  i  i  o
24  !*****************************************************************************
25  !                                                                            *
26  !  Calculation of the surface resistance according to the procedure given    *
27  !  in:                                                                       *
28  !  Wesely (1989): Parameterization of surface resistances to gaseous         *
29  !  dry deposition in regional-scale numerical models.                        *
30  !  Atmos. Environ. 23, 1293-1304.                                            *
31  !                                                                            *
32  !                                                                            *
33  !      AUTHOR: Andreas Stohl, 19 May 1995                                    *
34  !                                                                            *
35  !*****************************************************************************
36  !                                                                            *
37  ! Variables:                                                                 *
38  !                                                                            *
39  ! reldiff(maxspec)  diffusivity of H2O/diffusivity of component i            *
40  ! gr [W/m2]       global radiation                                           *
41  ! i               index of seasonal category                                 *
42  ! j               index of landuse class                                     *
43  ! ldep(maxspec)          1, if deposition shall be calculated for species i  *
44  ! nc                   actual number of chemical components                  *
45  ! rcl(maxspec,5,8) [s/m] Lower canopy resistance                             *
46  ! rgs(maxspec,5,8) [s/m] Ground resistance                                   *
47  ! rlu(maxspec,5,8) [s/m] Leaf cuticular resistance                           *
48  ! rm(maxspec) [s/m]      Mesophyll resistance                                *
49  ! t [C]           temperature                                                *
50  !                                                                            *
51  !*****************************************************************************
52
53  use par_mod
54  use com_mod
55
56  implicit none
57
58  integer :: i,j,ic,nc
59  real :: gr,rh,rr,t,rs,rsm,corr,rluc,rclc,rgsc,rdc,rluo
60  real :: rc(maxspec)
61
62
63  ! Compute stomatal resistance
64  !****************************
65  ! Sabine Eckhardt, Dec 06: use 1E25 instead of 99999. for infinite res.
66
67  if ((t.gt.0.).and.(t.lt.40.)) then
68    rs=ri(i,j)*(1.+(200./(gr+0.1))**2)*(400./(t*(40.-t)))
69  else
70    rs=1.E25
71  !  rs=99999.
72  endif
73
74
75  ! Correct stomatal resistance for effect of dew and rain
76  !*******************************************************
77
78  if ((rh.gt.0.9).or.(rr.gt.0.)) rs=rs*3.
79
80  ! Compute the lower canopy resistance
81  !************************************
82
83  rdc=100.*(1.+1000./(gr+10.))
84
85
86  corr=1000.*exp(-1.*t-4.)
87  do ic=1,nc
88    if (reldiff(ic).gt.0.) then
89
90  ! Compute combined stomatal and mesophyll resistance
91  !***************************************************
92
93      rsm=rs*reldiff(ic)+rm(ic)
94
95  ! Correct leaf cuticular, lower canopy and ground resistance
96  !***********************************************************
97
98      rluc=rlu(ic,i,j)+corr
99      rclc=rcl(ic,i,j)+corr
100      rgsc=rgs(ic,i,j)+corr
101
102  ! Correct leaf cuticular resistance for effect of dew and rain
103  !*************************************************************
104
105      if (rr.gt.0.) then
106        rluo=1./(1./1000.+1./(3.*rluc))
107        rluc=1./(1./(3.*rluc)+1.e-7*henry(ic)+f0(ic)/rluo)
108      else if (rh.gt.0.9) then
109        rluo=1./(1./3000.+1./(3.*rluc))
110        rluc=1./(1./(3.*rluc)+1.e-7*henry(ic)+f0(ic)/rluo)
111      endif
112
113  ! Combine resistances to give total resistance
114  !*********************************************
115
116      rc(ic)=1./(1./rsm+1./rluc+1./(rdc+rclc)+1./(rac(i,j)+rgsc))
117  ! Sabine Eckhardt, Dec 06: avoid possible excessively high vdep
118      if (rc(ic).lt.10.) rc(ic)=10.
119    endif
120  end do
121
122end subroutine getrc
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG