source: branches/jerome/src_flexwrf_v3.1/hanna_short.f90 @ 16

Last change on this file since 16 was 16, checked in by jebri, 10 years ago

sources for flexwrf v3.1

File size: 4.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 hanna_short(z, &
23  ust,wst,ol,h,zeta,sigu,sigv,tlu,tlv,tlw, &
24  sigw,dsigwdz,dsigw2dz)
25
26  !                       i
27  !*****************************************************************************
28  !                                                                            *
29  !   Computation of \sigma_i and \tau_L based on the scheme of Hanna (1982)   *
30  !                                                                            *
31  !   Author: A. Stohl                                                         *
32  !                                                                            *
33  !   4 December 1997                                                          *
34  !                                                                            *
35  !*****************************************************************************
36  !                                                                            *
37  ! Variables:                                                                 *
38  ! dsigwdz [1/s]     vertical gradient of sigw                                *
39  ! ol [m]            Obukhov length                                           *
40  ! sigu, sigv, sigw  standard deviations of turbulent velocity fluctuations   *
41  ! tlu [s]           Lagrangian time scale for the along wind component.      *
42  ! tlv [s]           Lagrangian time scale for the cross wind component.      *
43  ! tlw [s]           Lagrangian time scale for the vertical wind component.   *
44  ! ust, ustar [m/s]  friction velocity                                        *
45  ! wst, wstar [m/s]  convective velocity scale                                *
46  !                                                                            *
47  !*****************************************************************************
48
49  use par_mod
50  use com_mod
51!  use hanna_mod
52
53  implicit none
54
55  real :: z
56  real :: ust,wst,ol,h,zeta,sigu,sigv,tlu,tlv,tlw
57  real :: sigw,dsigwdz,dsigw2dz
58
59
60
61  !**********************
62  ! 1. Neutral conditions
63  !**********************
64
65  if (h/abs(ol).lt.1.) then
66    ust=max(1.e-4,ust)
67    sigw=1.3*exp(-2.e-4*z/ust)
68    dsigwdz=-2.e-4*sigw
69    sigw=sigw*ust+1.e-2
70    tlw=0.5*z/sigw/(1.+1.5e-3*z/ust)
71
72
73  !***********************
74  ! 2. Unstable conditions
75  !***********************
76
77  else if (ol.lt.0.) then
78
79
80  ! Determine sigmas
81  !*****************
82
83    sigw=sqrt(1.2*wst**2*(1.-.9*zeta)*zeta**0.66666+ &
84         (1.8-1.4*zeta)*ust**2)+1.e-2
85    dsigwdz=0.5/sigw/h*(-1.4*ust**2+wst**2* &
86         (0.8*max(zeta,1.e-3)**(-.33333)-1.8*zeta**0.66666))
87
88
89  ! Determine average Lagrangian time scale
90  !****************************************
91
92    if (z.lt.abs(ol)) then
93      tlw=0.1*z/(sigw*(0.55-0.38*abs(z/ol)))
94    else if (zeta.lt.0.1) then
95      tlw=0.59*z/sigw
96    else
97      tlw=0.15*h/sigw*(1.-exp(-5*zeta))
98    endif
99
100
101  !*********************
102  ! 3. Stable conditions
103  !*********************
104
105  else
106    sigw=1.e-2+1.3*ust*(1.-zeta)
107    dsigwdz=-1.3*ust/h
108    tlw=0.1*h/sigw*zeta**0.8
109  endif
110
111
112  tlu=max(10.,tlu)
113  tlv=max(10.,tlv)
114  tlw=max(30.,tlw)
115  if (dsigwdz.eq.0.) dsigwdz=1.e-10
116
117end subroutine hanna_short
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG