source: flexpart.git/src/getrb.f90

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bug
Last change on this file was 92fab65, checked in by Ignacio Pisso <ip@…>, 4 years ago

add SPDX-License-Identifier to all .f90 files

  • Property mode set to 100644
File size: 2.1 KB
Line 
1! SPDX-FileCopyrightText: FLEXPART 1998-2019, see flexpart_license.txt
2! SPDX-License-Identifier: GPL-3.0-or-later
3
4subroutine getrb(nc,ustar,nyl,diffh2o,reldiff,rb)
5  !                 i    i    i     i       i    o
6  !*****************************************************************************
7  !                                                                            *
8  !  Calculation of the quasilaminar sublayer resistance to dry deposition.    *
9  !                                                                            *
10  !      AUTHOR: Andreas Stohl, 20 May 1995                                    *
11  !                                                                            *
12  !*****************************************************************************
13  !                                                                            *
14  ! Variables:                                                                 *
15  ! rb(ncmax)       sublayer resistance                                        *
16  ! schmidt         Schmidt number                                             *
17  ! ustar [m/s]     friction velocity                                          *
18  ! diffh20 [m2/s]  diffusivity of water vapor in air                          *
19  ! reldiff         diffusivity relative to H2O                                *
20  !                                                                            *
21  ! Constants:                                                                 *
22  ! karman          von Karman constant                                        *
23  ! pr              Prandtl number                                             *
24  !                                                                            *
25  !*****************************************************************************
26
27  use par_mod
28
29  implicit none
30
31  real :: ustar,diffh2o,rb(maxspec),schmidt,nyl
32  real :: reldiff(maxspec)
33  integer :: ic,nc
34  real,parameter :: pr=0.72
35
36  do ic=1,nc
37    if (reldiff(ic).gt.0.) then
38      schmidt=nyl/diffh2o*reldiff(ic)
39      rb(ic)=2.0*(schmidt/pr)**0.67/(karman*ustar)
40    endif
41  end do
42
43end subroutine getrb
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG