The algebraic k-equation (Source File: tkealgebraic.F90)


INTERFACE:

    subroutine tkealgebraic(nlev,u_taus,u_taub,NN,SS)
DESCRIPTION:

This subroutine computes the turbulent kinetic energy based on (150), but using the local equilibrium assumption

$\displaystyle P+G-\epsilon=0 \quad .$ (169)

This statement can be re-expressed in the form

$\displaystyle k= (c_\mu^0)^{-3} \, l^2 ( c_\mu M^2 - c'_\mu N^2 ) \; , \quad$ (170)

were we used the expressions in (152) together with (43) and (44). The rate of dissipaton, $ \epsilon$, has been expressed in terms of $ l$ via (153). This equation has been implemented to update $ k$ in a diagnostic way. It is possible to compute the value of $ k$ as the weighted average of (170) and the value of $ k$ at the old timestep. The weighting factor is defined by the parameter c_filt. It is recommended to take this factor small (e.g. c_filt = 0.2) in order to reduce the strong oscillations associated with this scheme, and to couple it with an algebraically prescribed length scale with the length scale limitation active (length_lim=.true. in gotmturb.nml, see Galperin et al. (1988)).


USES:

    use turbulence,   only: tke,tkeo,L,k_min
    use turbulence,   only: cmue2,cde,cmue1,cm0
 
    IMPLICIT NONE
INPUT PARAMETERS:
 
    number of vertical layers
    integer,  intent(in)                :: nlev
 
    surface and bottom
    friction velocity (m/s)
    REALTYPE, intent(in)                :: u_taus,u_taub
 
    square of shear and buoyancy
    frequency (1/s^2)
    REALTYPE, intent(in)                :: NN(0:nlev),SS(0:nlev)
DEFINED PARAMETERS:
    REALTYPE , parameter                :: c_filt=1.0
REVISION HISTORY:
    Original author(s): Hans Burchard & Karsten Bolding
    $Log: tkealgebraic.F90,v $
    Revision 1.7  2007-01-06 11:49:15  kbk
    namelist file extension changed .inp --> .nml
    Revision 1.6  2005/06/27 13:44:07  kbk
    modified + removed traling blanks
    Revision 1.5  2003/03/28 09:20:35  kbk
    added new copyright to files
    Revision 1.4  2003/03/28 08:29:16  kbk
    removed tabs
    Revision 1.3  2003/03/10 09:02:05  gotm
    Added new Generic Turbulence Model + improved documentation and cleaned up code
    Revision 1.2  2002/02/08 08:59:58  gotm
    Revision 1.1.1.1  2001/02/12 15:55:58  gotm
    initial import into CVS

Karsten Bolding 2012-01-24