source: trunk/src/hanna1.f90 @ 28

Last change on this file since 28 was 4, checked in by mlanger, 11 years ago
File size: 5.1 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 hanna1(z)
23  !                  i
24  !*****************************************************************************
25  !                                                                            *
26  !   Computation of \sigma_i and \tau_L based on the scheme of Hanna (1982)   *
27  !                                                                            *
28  !   Author: A. Stohl                                                         *
29  !                                                                            *
30  !   4 December 1997                                                          *
31  !                                                                            *
32  !*****************************************************************************
33  !                                                                            *
34  ! Variables:                                                                 *
35  ! dsigwdz [1/s]     vertical gradient of sigw                                *
36  ! ol [m]            Obukhov length                                           *
37  ! sigu, sigv, sigw  standard deviations of turbulent velocity fluctuations   *
38  ! tlu [s]           Lagrangian time scale for the along wind component.      *
39  ! tlv [s]           Lagrangian time scale for the cross wind component.      *
40  ! tlw [s]           Lagrangian time scale for the vertical wind component.   *
41  ! ust, ustar [m/s]  friction velocity                                        *
42  ! wst, wstar [m/s]  convective velocity scale                                *
43  !                                                                            *
44  !*****************************************************************************
45
46  use par_mod
47  use com_mod
48  use hanna_mod
49
50  implicit none
51
52  real :: z,s1,s2
53
54
55
56  !**********************
57  ! 1. Neutral conditions
58  !**********************
59
60  if (h/abs(ol).lt.1.) then
61
62    ust=max(1.e-4,ust)
63    sigu=2.0*ust*exp(-3.e-4*z/ust)
64    sigu=max(sigu,1.e-5)
65    sigv=1.3*ust*exp(-2.e-4*z/ust)
66    sigv=max(sigv,1.e-5)
67    sigw=sigv
68    dsigw2dz=-6.76e-4*ust*exp(-4.e-4*z/ust)
69    tlu=0.5*z/sigw/(1.+1.5e-3*z/ust)
70    tlv=tlu
71    tlw=tlu
72
73
74  !***********************
75  ! 2. Unstable conditions
76  !***********************
77
78  else if (ol.lt.0.) then
79
80
81  ! Determine sigmas
82  !*****************
83
84    sigu=ust*(12-0.5*h/ol)**0.33333
85    sigu=max(sigu,1.e-6)
86    sigv=sigu
87
88    if (zeta.lt.0.03) then
89      sigw=0.96*wst*(3*zeta-ol/h)**0.33333
90      dsigw2dz=1.8432*wst*wst/h*(3*zeta-ol/h)**(-0.33333)
91    else if (zeta.lt.0.4) then
92      s1=0.96*(3*zeta-ol/h)**0.33333
93      s2=0.763*zeta**0.175
94      if (s1.lt.s2) then
95        sigw=wst*s1
96        dsigw2dz=1.8432*wst*wst/h*(3*zeta-ol/h)**(-0.33333)
97      else
98        sigw=wst*s2
99        dsigw2dz=0.203759*wst*wst/h*zeta**(-0.65)
100      endif
101    else if (zeta.lt.0.96) then
102      sigw=0.722*wst*(1-zeta)**0.207
103      dsigw2dz=-.215812*wst*wst/h*(1-zeta)**(-0.586)
104    else if (zeta.lt.1.00) then
105      sigw=0.37*wst
106      dsigw2dz=0.
107    endif
108    sigw=max(sigw,1.e-6)
109
110
111  ! Determine average Lagrangian time scale
112  !****************************************
113
114    tlu=0.15*h/sigu
115    tlv=tlu
116    if (z.lt.abs(ol)) then
117      tlw=0.1*z/(sigw*(0.55-0.38*abs(z/ol)))
118    else if (zeta.lt.0.1) then
119      tlw=0.59*z/sigw
120    else
121      tlw=0.15*h/sigw*(1.-exp(-5*zeta))
122    endif
123
124
125  !*********************
126  ! 3. Stable conditions
127  !*********************
128
129  else
130    sigu=2.*ust*(1.-zeta)
131    sigv=1.3*ust*(1.-zeta)
132    sigu=max(sigu,1.e-6)
133    sigv=max(sigv,1.e-6)
134    sigw=sigv
135    dsigw2dz=3.38*ust*ust*(zeta-1.)/h
136    tlu=0.15*h/sigu*(sqrt(zeta))
137    tlv=0.467*tlu
138    tlw=0.1*h/sigw*zeta**0.8
139  endif
140
141
142
143
144  tlu=max(10.,tlu)
145  tlv=max(10.,tlv)
146  tlw=max(30.,tlw)
147
148
149end subroutine hanna1
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG