TROENKz Subroutine

private subroutine TROENKz(z, h, ustar, invL, Kdef, Kz, phih)

Arguments

Type IntentOptional Attributes Name
real, intent(in) :: z
real, intent(in) :: h
real, intent(in) :: ustar
real, intent(in) :: invL
real, intent(in) :: Kdef
real, intent(out) :: Kz
real, intent(out) :: phih

Called by

proc~~troenkz~~CalledByGraph proc~troenkz TROENKz proc~uemep_set_dispersion_sigma_kz uEMEP_set_dispersion_sigma_Kz proc~uemep_set_dispersion_sigma_kz->proc~troenkz proc~uemep_subgrid_deposition uEMEP_subgrid_deposition proc~uemep_subgrid_deposition->proc~uemep_set_dispersion_sigma_kz proc~uemep_subgrid_dispersion uEMEP_subgrid_dispersion proc~uemep_subgrid_dispersion->proc~uemep_set_dispersion_sigma_kz proc~uemep_subgrid_dispersion_integral uEMEP_subgrid_dispersion_integral proc~uemep_subgrid_dispersion->proc~uemep_subgrid_dispersion_integral proc~uemep_subgrid_dispersion_integral->proc~uemep_set_dispersion_sigma_kz program~uemep uEMEP program~uemep->proc~uemep_subgrid_deposition program~uemep->proc~uemep_subgrid_dispersion

Source Code

    subroutine TROENKz(z,h,ustar,invL,Kdef,Kz,phih)
        !Vertical dispersion routine as described in:
        !Troen, I., & Mahrt, L. (1986). A Simple Model of the Atmospheric Boundary
        !Layer:
        !Sensitivity to Surface Evaporation. Boundary-Layer Meteorology, 37, 129-148.
        !https://doi.org/10.1007/BF00122760

        implicit none

        real, intent(in) :: z    ! height
        real, intent(in) :: h    ! Boundary layer depth
        real, intent(in) :: ustar, invL, Kdef !  u*, 1/L, default Kz
        real, intent(out) :: Kz,phih
        real :: ws
        real :: kappa=0.4
        real :: zsurf

        !Height of the surface layer. Stability function for unstable calculated
        !here for z>zsurf
        zsurf=0.1*h

        if ( z < h ) then
            if ( invL < 0 ) then
                !phih=(1-7.*min(z,zsurf)*invL)**(-1./3.) !Original in
                !Troen and Mahrt for phim, so no Prandtl number
                !phih=(1-16.*min(z,zsurf)*invL)**(-1./2.) !As in Garratt and Obrien for phih, so with Prandtl number
                phih=(1-16.*z*invL)**(-1./2.) !As in Garratt and Obrien for phih, so with Prandtl number
            else
                phih=1+5.*z*invL !As in Garratt, Prandtl number is 1 in stable boundary layer
            endif
            ws=ustar/phih
            Kz=kappa*ws*z*(1.-z/h)**2

        else
            Kz =  Kdef
        endif

    end subroutine TROENKz