source: branches/jerome/src_flexwrf_v3.1/tke_partition_hanna.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: 7.5 KB
Line 
1!***********************************************************************
2!* Copyright 2012,2013                                                *
3!* Jerome Brioude, Delia Arnold, Andreas Stohl, Wayne Angevine,       *
4!* John Burkhart, Massimo Cassiani, Adam Dingwell, Richard C Easter, Sabine Eckhardt,*
5!* Stephanie Evan, Jerome D Fast, Don Morton, Ignacio Pisso,          *
6!* Petra Seibert, Gerard Wotawa, Caroline Forster, Harald Sodemann,   *
7!*                                                                     *
8!* This file is part of FLEXPART WRF                                   *
9!*                                                                     *
10!* FLEXPART is free software: you can redistribute it and/or modify    *
11!* it under the terms of the GNU General Public License as published by*
12!* the Free Software Foundation, either version 3 of the License, or   *
13!* (at your option) any later version.                                 *
14!*                                                                     *
15!* FLEXPART is distributed in the hope that it will be useful,         *
16!* but WITHOUT ANY WARRANTY; without even the implied warranty of      *
17!* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
18!* GNU General Public License for more details.                        *
19!*                                                                     *
20!* You should have received a copy of the GNU General Public License   *
21!* along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
22!***********************************************************************
23
24      subroutine tke_partition_hanna(z, &
25   ust,wst,ol,h,zeta,sigu,sigv,tlu,tlv,tlw, &
26   sigw,dsigwdz,dsigw2dz,uprof,vprof,tkeprof,pttprof,indz,indzp)
27!                      i
28!*******************************************************************************
29!                                                                              *
30!     Computation of \sigma_u,v,w,dsigwdz, and dsigw2dz based on TKE from WRF  *
31!                                                                              *
32!     Author: A. Stohl                                                         *
33!                                                                              *
34!     4 December 1997                                                          *
35!                                                                              *
36!*******************************************************************************
37!                                                                              *
38! Variables:                                                                   *
39! dsigwdz [1/s]     vertical gradient of sigw                                  *
40! ol [m]            Obukhov length                                             *
41! sigu, sigv, sigw  standard deviations of turbulent velocity fluctuations     *
42! tlu [s]           Lagrangian time scale for the along wind component.        *
43! tlv [s]           Lagrangian time scale for the cross wind component.        *
44! tlw [s]           Lagrangian time scale for the vertical wind component.     *
45! ust, ustar [m/s]  friction velocity                                          *
46! wst, wstar [m/s]  convective velocity scale                                  *
47!                                                                              *
48!*******************************************************************************
49! 12 JUNE 2007  compute sigu,sigv,sigw,dsigwdz, and dsigw2dz from TKE          *
50! fu2               fraction for u2
51! fv2               fraction for v2
52! fw2               fraction for w2
53!
54! 25 JUNE 2007 merged from hanna_tke and tke_partition.f
55
56!      include 'includepar'
57!      include 'includecom'
58!      include 'includehanna'
59!      include 'includeinterpol'
60  use par_mod
61  use com_mod
62!  use hanna_mod
63!  use interpol_mod
64   implicit none
65      real :: corr,z,zzz,hanna_sigu,hanna_sigv,hanna_sigw,dz,dz1,dz2
66      real :: fu2(2),fv2(2),fw2(2)
67      real :: siguprof(2),sigvprof(2),sigwprof(2)
68      integer :: k,indz,indzp
69     real :: uprof(nzmax),vprof(nzmax),tkeprof(nzmax),pttprof(nzmax)
70    real :: ust,wst,ol,h,zeta,sigu,sigv,tlu,tlv,tlw
71    real :: sigw,dsigwdz,dsigw2dz
72
73       do k=1,2
74         fu2(k)=0.333
75         fv2(k)=0.333
76         fw2(k)=0.333       
77         siguprof(k)=0.0
78         sigvprof(k)=0.0
79         sigwprof(k)=0.0
80       enddo
81
82! determine fraction
83       do k=1,2                                   ! k fraction     
84!**********************
85! 1. Neutral conditions
86!**********************
87        if (k .eq. 1) zzz=height(indz)
88        if (k .eq. 2) zzz=height(indzp)
89      if (h/abs(ol).lt.1.) then
90        ust=max(1.e-4,ust)
91        corr=height(indz+k-1)/ust
92        hanna_sigu=1.e-2+2.0*ust*exp(-3.e-4*corr)
93        hanna_sigw=1.e-2+1.3*ust*exp(-2.e-4*corr)
94        hanna_sigv=hanna_sigw
95
96!***********************
97! 2. Unstable conditions
98!***********************
99
100      else if (ol.lt.0.) then
101
102
103! Determine sigmas
104!*****************
105
106        hanna_sigu=1.e-2+ust*(12-0.5*h/ol)**0.33333
107        hanna_sigv=hanna_sigu
108        hanna_sigw=sqrt(1.2*wst**2*(1.-.9*zeta)*height(indz+k-1)/h &
109       **0.66666+(1.8-1.4*height(k)/h)*ust**2)+1.e-2
110
111!*********************
112! 3. Stable conditions
113!*********************
114      else
115        hanna_sigu=1.e-2+2.*ust*(1.-height(indz+k-1)/h)
116        hanna_sigv=1.e-2+1.3*ust*(1.-height(indz+k-1)/h)
117        hanna_sigw=hanna_sigv
118      endif
119        fu2(k)=hanna_sigu**2/(hanna_sigu**2+hanna_sigv**2+hanna_sigw**2)
120        fv2(k)=hanna_sigv**2/(hanna_sigu**2+hanna_sigv**2+hanna_sigw**2)
121        fw2(k)=hanna_sigw**2/(hanna_sigu**2+hanna_sigv**2+hanna_sigw**2)
122      enddo                                          !k fraction
123
124!- compute sigu,v,w
125      do k=1,2             ! siguprof
126        siguprof(k)=max(sqrt(2.0*tkeprof(indz+k-1)*fu2(k)),1.e-2)
127        sigvprof(k)=max(sqrt(2.0*tkeprof(indz+k-1)*fv2(k)),1.e-2)
128        sigwprof(k)=max(sqrt(2.0*tkeprof(indz+k-1)*fw2(k)),1.e-2)
129!C         write(*,*)'z=',height(indz+k-1),'tke=', tkeprof(indz+k-1)
130      enddo                 ! siguprof
131!         write(*,*)'tkeprof=',(tkeprof(k),k=1,nz)
132!- interpolate sigu,sigv, sigw
133          dz=1./(height(indzp)-height(indz))
134          dz1=(z - height(indz))*dz
135          dz2=(height(indzp)-z)*dz
136           
137          sigu=dz1*siguprof(2)+dz2*siguprof(1)
138          sigv=dz1*sigvprof(2)+dz2*sigvprof(1)
139          sigw=dz1*sigwprof(2)+dz2*sigwprof(1)
140         
141          dsigwdz=max(1.e-10,(sigwprof(2)-sigwprof(1))*dz )         
142          dsigw2dz=max(1.e-10,(sigwprof(2)**2-sigwprof(1)**2)*dz )
143
144!-- compute length scales based on hanna(1982)
145!*************************
146! 1. Neutral conditions
147!**********************
148 
149      if (h/abs(ol).lt.1.) then
150        tlu=0.5*z/sigw/(1.+1.5e-3*corr)
151        tlv=tlu
152        tlw=tlu
153 
154!***********************
155! 2. Unstable conditions
156!***********************
157 
158      else if (ol.lt.0.) then
159 
160! Determine average Lagrangian time scale
161!****************************************
162 
163        tlu=0.15*h/sigu
164        tlv=tlu
165        if (z.lt.abs(ol)) then
166          tlw=0.1*z/(sigw*(0.55-0.38*abs(z/ol)))
167        else if (zeta.lt.0.1) then
168          tlw=0.59*z/sigw
169        else
170          tlw=0.15*h/sigw*(1.-exp(-5*zeta))
171        endif
172 
173 
174!*********************
175! 3. Stable conditions
176!*********************
177 
178      else
179        tlu=0.15*h/sigu*(sqrt(zeta))
180        tlv=0.467*tlu
181        tlw=0.1*h/sigw*zeta**0.8
182      endif
183 
184      tlu=max(10.,tlu)
185      tlv=max(10.,tlv)
186      tlw=max(30.,tlw)
187 
188
189
190
191
192      end subroutine tke_partition_hanna
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG