module obl_k_epsilon
!< @brief k-epsilon scheme module
! --------------------------------------------------------------------------------
! modules used
! --------------------------------------------------------------------------------
#ifdef USE_CONFIG_PARSER
use iso_c_binding, only: C_NULL_CHAR
use config_parser
#endif
use obl_math
use obl_grid
use obl_bc
use obl_turb_common
! directives list
! --------------------------------------------------------------------------------
implicit none
private
! public interface
! --------------------------------------------------------------------------------
public :: init_turbulence_closure
public :: advance_turbulence_closure
public :: define_stability_functions
public :: get_eddy_viscosity, get_eddy_diffusivity
public :: get_eddy_viscosity_vec, get_eddy_diffusivity_vec
public :: TKE_init, EPS_init
public :: TKE_bc, EPS_bc
public :: set_config_param
! --------------------------------------------------------------------------------
!> @brief k-epsilon parameters
type, public :: kepsilonParamType
real :: C_mu = 0.09 !< momentum stability function constant
real :: C1_eps = 1.44 !< eps-eq.: production constant
real :: C2_eps = 1.92 !< eps-eq.: dissipation constant
real :: C3_eps = -0.40 !< eps-eq.: buoyancy constant (-0.4 stable stratification or 1.14 unstable stratification)
real :: ScT_TKE = 1.0 !< Schmidt number for TKE
real :: ScT_EPS = 1.11 !< Scmidt number for eps
real :: TKE_min = 0.000001 !< minimal valur for TKE, [J/kg]
real :: EPS_min = (0.000001)**(3.0/2.0) !< minimal value for epsilon (TKE dissipation rate), [W/kg]
real :: kappa = 0.4 !< von Karman constant
real :: PrT0 = 1.25 !< turbulent Prandtl number
real :: c_S2_min = 1e-5 !< min shear frequency, [(1/s)**2]
end type
contains
! --------------------------------------------------------------------------------
subroutine define_stability_functions(param, bc, grid)
!< @brief advance stability functions: N**2, S**2, Ri-gradient
! ----------------------------------------------------------------------------
use obl_state
type(kepsilonParamType), intent(in) :: param
type(oblBcType), intent(in) :: bc
type (gridDataType), intent(in) :: grid
! --------------------------------------------------------------------------------
call get_N2(N2, Rho, bc%Rho_dyn0, bc%Rho_dynH, &
param%kappa, param%PrT0, grid)
call get_S2(S2, U, V, bc%U_dyn0, bc%U_dynH, &
param%kappa, grid)
call get_Ri_gradient(Ri_grad, N2, S2, param%c_S2_min, grid)
end subroutine
! --------------------------------------------------------------------------------
subroutine init_turbulence_closure(param, bc, grid)
!< @brief advance turbulence closure
! ----------------------------------------------------------------------------
use obl_state
type(kepsilonParamType), intent(in) :: param
type(oblBcType), intent(in) :: bc
type (gridDataType), intent(in) :: grid
! --------------------------------------------------------------------------------
call TKE_init(TKE, param, grid%cz)
call EPS_init(EPS, param, grid%cz, grid%height)
call TKE_bc(TKE, &
bc%U_dyn0, bc%U_dynH, param, grid%cz)
call EPS_bc(EPS, &
bc%U_dyn0, bc%U_dynH, param, grid%dz, grid%cz)
call get_eddy_viscosity(Km, TKE, EPS, param, grid)
call get_eddy_diffusivity(Kh, Km, param, grid)
end subroutine
! --------------------------------------------------------------------------------
subroutine advance_turbulence_closure(param, bc, grid, dt)
!< @brief advance turbulence closure
! ----------------------------------------------------------------------------
use obl_state
type(kepsilonParamType), intent(in) :: param
type(oblBcType), intent(in) :: bc
type (gridDataType), intent(in) :: grid
real, intent(in) :: dt
! --------------------------------------------------------------------------------
call get_TKE_production(TKE_production, Km, S2, grid)
call get_TKE_buoyancy(TKE_buoyancy, Kh, N2, grid)
call TKE_bc(TKE, &
bc%U_dyn0, bc%U_dynH, param, grid%cz)
call EPS_bc(EPS, &
bc%U_dyn0, bc%U_dynH, param, grid%dz, grid%cz)
call get_eddy_viscosity(Km, TKE, EPS, param, grid)
call get_eddy_diffusivity(Kh, Km, param, grid)
call TKE_bc(TKE, &
bc%U_dyn0, bc%U_dynH, param, grid%cz)
call eps_bc(EPS, &
bc%U_dyn0, bc%U_dynH, param, grid%dz, grid%cz)
!call solve_TKE_eq_semiimplicit (TKE, Km, grid%cz, grid%dz, dt, TKE_production, TKE_buoyancy, eps, param)
call solve_TKE_eq (TKE, Km, grid%cz, grid%dz, dt, TKE_production, TKE_buoyancy, EPS, param)
call limit_min_array(TKE, param%TKE_min, grid%cz)
call solve_eps_eq_semiimplicit (EPS, Km, grid%cz, grid%dz, dt, TKE_production, TKE_buoyancy, TKE, param)
call limit_min_array(EPS, param%eps_min, grid%cz)
end subroutine
subroutine TKE_init(TKE, param, nz)
!< @brief TKE initialization
type(kepsilonParamType), intent(in) :: param
integer, intent(in) :: nz !< number of z-steps
real, intent(out) :: TKE(nz) !< TKE, [J/kg]
integer :: k !< counter
do k = 1, nz
TKE(k) = 1.0 / (param%C_mu)**(1.0/2.0) * 0.001 * 0.001
end do
end subroutine
subroutine EPS_init(eps, param, nz, z)
!< @brief TKE dissipation rate initialization
type(kepsilonParamType), intent(in) :: param
integer, intent(in) :: nz !< number of z-steps
real, intent(in) :: z !< depth, [m]
real, intent(out) :: eps(nz) !< TKE dissipation rate, [W/kg]
integer :: k !< counter
do k = 1, nz
eps(k) = 0.001 * 0.001 * 0.001 / (param%kappa * z) !not full z, but current z - must be fixed
end do
end subroutine
! ----------------------------------------------------------------------------
subroutine TKE_bc(TKE, U_dyn0, U_dynH, param, cz)
!< @brief TKE boundary conditions
! ----------------------------------------------------------------------------
type(kepsilonParamType), intent(in) :: param
real, intent(in) :: U_dyn0, U_dynH !< dynamic velocity, [(m/s)]
integer, intent(in) :: cz !< grid size
real, intent(out) :: TKE(cz) !< TKE, [(m/s)**2]
! ----------------------------------------------------------------------------
TKE(1) = (1.0 / sqrt(param%C_mu)) * U_dyn0 * U_dyn0
TKE(cz) = (1.0 / sqrt(param%C_mu)) * U_dynH * U_dynH
end subroutine TKE_bc
! ----------------------------------------------------------------------------
subroutine EPS_bc(EPS, U_dyn0, U_dynH, param, dz, cz)
!< @brief epsilon boundary conditions
! ----------------------------------------------------------------------------
type(kepsilonParamType), intent(in) :: param
real, intent(in) :: U_dyn0, U_dynH !< dynamic velocity, [(m/s)]
integer, intent(in) :: cz !< grid size
real, intent(in) :: dz(cz) !< grid steps
real, intent(out) :: EPS(cz) !< dissipation rate, [m^2/s^3]
EPS(1) = (U_dyn0 * U_dyn0 * U_dyn0) / (param%kappa * 0.5 * dz(1))
EPS(cz) = (U_dynH * U_dynH * U_dynH) / (param%kappa * 0.5 * dz(cz))
end subroutine EPS_bc
! --------------------------------------------------------------------------------
subroutine get_eddy_viscosity(Km, TKE, EPS, param, grid)
!< @brief calculate eddy viscosity
! ----------------------------------------------------------------------------
type(kepsilonParamType), intent(in) :: param
type (gridDataType), intent(in) :: grid
real, dimension(grid%cz), intent(out) :: Km !< eddy viscosity, [m**2 / s]
real, dimension(grid%cz), intent(in) :: TKE, EPS !<
! ----------------------------------------------------------------------------
call get_eddy_viscosity_vec(Km, TKE, EPS, param, grid%cz)
end subroutine get_eddy_viscosity
subroutine get_eddy_viscosity_vec(Km, TKE, EPS, param, cz)
!< @brief calculate eddy viscosity
! ----------------------------------------------------------------------------
type(kepsilonParamType), intent(in) :: param
integer, intent(in) :: cz !< grid size
real, dimension(cz), intent(out) :: Km !< eddy viscosity, [m**2 / s]
real, dimension(cz), intent(in) :: TKE, EPS !<
integer :: k
! ----------------------------------------------------------------------------
do k = 1, cz
if (eps(k) < param%eps_min) then
Km(k) = param%C_mu * TKE(k) * TKE(k) / param%eps_min
else
Km(k) = param%C_mu * TKE(k) * TKE(k) / eps(k)
endif
end do
end subroutine get_eddy_viscosity_vec
! --------------------------------------------------------------------------------
subroutine get_eddy_diffusivity(Kh, Km, param, grid)
!< @brief calculate eddy diffusivity
! ----------------------------------------------------------------------------
type(kepsilonParamType), intent(in) :: param
type (gridDataType), intent(in) :: grid
real, dimension(grid%cz), intent(out) :: Kh !< eddy diffusivity, [m**2 / s]
real, dimension(grid%cz), intent(in) :: Km !< eddy viscosity, [m**2 / s]
! ----------------------------------------------------------------------------
call get_eddy_diffusivity_vec(Kh, Km, param, grid%cz)
end subroutine get_eddy_diffusivity
subroutine get_eddy_diffusivity_vec(Kh, Km, param, cz)
!< @brief calculate eddy diffusivity
! ----------------------------------------------------------------------------
type(kepsilonParamType), intent(in) :: param
integer, intent(in) :: cz !< grid size
real, dimension(cz), intent(out) :: Kh !< eddy diffusivity, [m**2 / s]
real, dimension(cz), intent(in) :: Km !< eddy viscosity, [m**2 / s]
integer :: k
! ----------------------------------------------------------------------------
do k = 1, cz
Kh(k) = Km(k) / param%PrT0
end do
end subroutine get_eddy_diffusivity_vec
! --------------------------------------------------------------------------------
subroutine solve_TKE_eq (Field, Km, nz, dz, dt, P_TKE, B_TKE, eps, param)
!< @brief solver for TKE equation (explicit)
type(kepsilonParamType), intent(in) :: param
integer, intent (in) :: nz !< number of steps in z (depth), [m]
real, intent (in) :: dt !< time step [s]
real, intent (in) :: dz(nz) !< depth(z) step [m]
real, intent(in) :: Km(nz)
real :: Kvisc(nz) !< turbulent transfer coefficient [m^2/s]
real :: G(nz) !< explicit part
real :: a(nz), b(nz), c(nz), f(nz) !< coefficients for tridiagonal matrix
real :: Field(nz) !< value of unknown field in the point
integer :: k !< counter
real :: Kvisc_plus(nz), Kvisc_minus(nz) !<half-sums for turbulent transfer coefficient Kvisc
real :: dz_plus(nz), dz_minus(nz) !<half-sums for depth-step dz
real :: f_right(nz) !< right-hand side
real :: gamma !< variable determining the explicitness/implicitness of a scheme (from 0.0 to 1.0)
real :: P_TKE(nz) !< shear production of TKE, [m**2 / s**3]
real :: B_TKE(nz) !< buoyancy production of TKE, [m**2 / s**3]
real :: eps(nz) !< TKE dissipation rate, [W/kg]
Kvisc(:) = Km(:) / param%ScT_TKE
gamma = 0.0
do k = 1, nz
f_right(k) = P_TKE(k) + B_TKE(k) - eps(k)
end do
do k = 1, nz-1
Kvisc_plus(k) = (Kvisc(k)+Kvisc(k+1))/2.0
dz_plus(k) = (dz(k)+dz(k+1))/2.0
end do
do k = 2, nz
Kvisc_minus(k) = (Kvisc(k)+Kvisc(k-1))/2.0
dz_minus(k) = (dz(k)+dz(k-1))/2.0
end do
do k = 2, nz-1
G(k) = Kvisc_plus(k)/dz_plus(k) * (Field(k+1) - Field(k)) &
- Kvisc_minus(k)/dz_minus(k) * (Field(k) - Field(k-1))
end do
G(1) = Kvisc_plus(1)/dz_plus(1) * (Field(2) - Field(1))
G(nz) = - Kvisc_minus(nz)/dz_minus(nz) * (Field(nz) - Field(nz-1))
do k = 2, nz-1
a(k) = dt / dz(k) * Kvisc_minus(k)/dz_minus(k) * (1.0 - gamma)
b(k) = -1.0 - dt / dz(k) * Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma) &
- dt / dz(k) * Kvisc_minus(k)/dz_minus(k) * (1.0 - gamma)
c(k) = dt / dz(k) * Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma)
f(k) = - Field(k) - f_right(k) * dt &
- dt / dz(k) * G(k) * gamma
if (k == 2) then
f(k) = f(k) - a(k) * Field(k-1)
endif
if (k == nz-1) then
f(k) = f(k) - c(k) * Field(k+1)
endif
if (k == 2) a(k) = 0.0
if (k == nz-1) c(k) = 0.0
end do
call tma (Field(2:nz-1), a(2:nz-1), b(2:nz-1), c(2:nz-1), f(2:nz-1), nz - 2)
end subroutine solve_TKE_eq
subroutine solve_eps_eq (Field, Km, nz, dz, dt, P_TKE, B_TKE, TKE, param)
!< @brief solver for epsilon equation (explicit)
type(kepsilonParamType), intent(in) :: param
integer, intent (in) :: nz !< number of steps in z (depth), [m]
real, intent (in) :: dt !< time step [s]
real, intent (in) :: dz(nz) !< depth(z) step [m]
real, intent(in) :: Km(nz)
real :: Kvisc(nz) !< turbulent transfer coefficient [m^2/s]
real :: G(nz) !< explicit part
real :: a(nz), b(nz), c(nz), f(nz) !< coefficients for tridiagonal matrix
real :: Field(nz) !< value of unknown field in the point
integer :: k !< counter
real :: Kvisc_plus(nz), Kvisc_minus(nz) !<half-sums for turbulent transfer coefficient Kvisc
real :: dz_plus(nz), dz_minus(nz) !<half-sums for depth-step dz
real :: f_right(nz) !< right-hand side
real :: gamma !< variable determining the explicitness/implicitness of a scheme (from 0.0 to 1.0)
real :: P_TKE(nz) !< shear production of TKE, [m**2 / s**3]
real :: B_TKE(nz) !< buoyancy production of TKE, [m**2 / s**3]
real :: eps(nz) !< TKE dissipation rate, [W/kg]
real :: TKE(nz) !< TKE, [J/kg]
Kvisc(:) = Km(:) / param%ScT_EPS
gamma = 0.0
do k = 1, nz
f_right(k) = Field(k) / TKE(k) * (param%C1_eps * P_TKE(k) + param%C3_eps * B_TKE(k) - param%C2_eps * Field(k))
end do
do k = 1, nz-1
Kvisc_plus(k) = (Kvisc(k)+Kvisc(k+1))/2.0
dz_plus(k) = (dz(k)+dz(k+1))/2.0
end do
do k = 2, nz
Kvisc_minus(k) = (Kvisc(k)+Kvisc(k-1))/2.0
dz_minus(k) = (dz(k)+dz(k-1))/2.0
end do
do k = 2, nz-1
G(k) = Kvisc_plus(k)/dz_plus(k) * (Field(k+1) - Field(k)) &
- Kvisc_minus(K)/dz_minus(k) * (Field(k) - Field(k-1))
end do
G(1) = Kvisc_plus(1)/dz_plus(1) * (Field(2) - Field(1))
G(nz) = - Kvisc_minus(nz)/dz_minus(nz) * (Field(nz) - Field(nz-1))
do k = 2, nz-1
a(k) = dt / dz(k) * Kvisc_minus(k)/dz_minus(k) * (1.0 - gamma)
b(k) = -1.0 - dt / dz(k) * Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma) &
- dt / dz(k) * Kvisc_minus(k)/dz_minus(k) * (1.0 - gamma)
c(k) = dt / dz(k) * Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma)
f(k) = - Field(k) - f_right(k) * dt &
- dt / dz(k) * G(k) * gamma
if (k == 2) then
f(k) = f(k) - a(k) * Field(k-1)
endif
if (k == nz-1) then
f(k) = f(k) - c(k) * Field(k+1)
endif
if (k == 2) a(k) = 0.0
if (k == nz-1) c(k) = 0.0
end do
call tma (Field(2:nz-1), a(2:nz-1), b(2:nz-1), c(2:nz-1), f(2:nz-1), nz - 2)
end subroutine solve_eps_eq
subroutine solve_TKE_eq_semiimplicit (Field, Km, nz, dz, dt, P_TKE, B_TKE, eps, param)
!< @brief solver for TKE equation (semimplicit)
type(kepsilonParamType), intent(in) :: param
integer, intent (in) :: nz !< number of steps in z (depth), [m]
real, intent (in) :: dt !< time step [s]
real, intent (in) :: dz(nz) !< depth(z) step [m]
real, intent(in) :: Km(nz)
real :: Kvisc(nz) !< turbulent transfer coefficient [m^2/s]
real :: G(nz) !< explicit part
real :: a(nz), b(nz), c(nz), f(nz) !< coefficients for tridiagonal matrix
real :: Field(nz) !< value of unknown field in the point
integer :: k !< counter
real :: Kvisc_plus(nz), Kvisc_minus(nz) !<half-sums for turbulent transfer coefficient Kvisc
real :: dz_plus(nz), dz_minus(nz) !<half-sums for depth-step dz
real :: f_right(nz) !< right-hand side
real :: gamma !< variable determining the explicitness/implicitness of a scheme (from 0.0 to 1.0)
real :: P_TKE(nz) !< shear production of TKE, [m**2 / s**3]
real :: B_TKE(nz) !< buoyancy production of TKE, [m**2 / s**3]
real :: eps(nz) !< TKE dissipation rate, [W/kg]
Kvisc(:) = Km(:) / param%ScT_TKE
gamma = 0.0
do k = 1, nz
f_right(k) = P_TKE(k) + B_TKE(k)
end do
do k = 1, nz-1
Kvisc_plus(k) = (Kvisc(k)+Kvisc(k+1))/2.0
dz_plus(k) = (dz(k)+dz(k+1))/2.0
end do
do k = 2, nz
Kvisc_minus(k) = (Kvisc(k)+Kvisc(k-1))/2.0
dz_minus(k) = (dz(k)+dz(k-1))/2.0
end do
do k = 2, nz-1
G(k) = Kvisc_plus(k)/dz_plus(k) * (Field(k+1) - Field(k)) &
- Kvisc_minus(k)/dz_minus(k) * (Field(k) - Field(k-1))
end do
G(1) = Kvisc_plus(1)/dz_plus(1) * (Field(2) - Field(1))
G(nz) = - Kvisc_minus(nz)/dz_minus(nz) * (Field(nz) - Field(nz-1))
do k = 2, nz-1
a(k) = dt / dz(k) * Kvisc_minus(k)/dz_minus(k) * (1.0 - gamma)
b(k) = -1.0 - dt * eps(k)/Field(k) - dt / dz(k) * &
Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma) &
- dt / dz(k) * Kvisc_minus(k)/dz_minus(k) &
* (1.0 - gamma)
c(k) = dt / dz(k) * Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma)
f(k) = - Field(k) - f_right(k) * dt &
- dt / dz(k) * G(k) * gamma
if (k == 2) then
f(k) = f(k) - a(k) * Field(k-1)
endif
if (k == nz-1) then
f(k) = f(k) - c(k) * Field(k+1)
endif
if (k == 2) a(k) = 0.0
if (k == nz-1) c(k) = 0.0
end do
call tma (Field(2:nz-1), a(2:nz-1), b(2:nz-1), c(2:nz-1), f(2:nz-1), nz - 2)
end subroutine solve_TKE_eq_semiimplicit
subroutine solve_eps_eq_semiimplicit (Field, Km, nz, dz, dt, P_TKE, B_TKE, TKE, param)
!< @brief solver for epsilon equation (semiimplict)
type(kepsilonParamType), intent(in) :: param
integer, intent (in) :: nz !< number of steps in z (depth), [m]
real, intent (in) :: dt !< time step [s]
real, intent (in) :: dz(nz) !< depth(z) step [m]
real, intent(in) :: Km(nz)
real :: Kvisc(nz) !< turbulent transfer coefficient [m^2/s]
real :: G(nz) !< explicit part
real :: a(nz), b(nz), c(nz), f(nz) !< coefficients for tridiagonal matrix
real :: Field(nz) !< value of unknown field in the point
integer :: k !< counter
real :: Kvisc_plus(nz), Kvisc_minus(nz) !<half-sums for turbulent transfer coefficient Kvisc
real :: dz_plus(nz), dz_minus(nz) !<half-sums for depth-step dz
real :: f_right(nz) !< right-hand side
real :: gamma !< variable determining the explicitness/implicitness of a scheme (from 0.0 to 1.0)
real :: P_TKE(nz) !< shear production of TKE, [m**2 / s**3]
real :: B_TKE(nz) !< buoyancy production of TKE, [m**2 / s**3]
real :: eps(nz) !< TKE dissipation rate, [W/kg]
real :: TKE(nz) !< TKE, [J/kg]
Kvisc(:) = Km(:) / param%ScT_EPS
gamma = 0.0
do k = 1, nz
f_right(k) = Field(k) / TKE(k) * (param%C1_eps * P_TKE(k) + param%C3_eps * B_TKE(k))
end do
do k = 1, nz-1
Kvisc_plus(k) = (Kvisc(k)+Kvisc(k+1))/2.0
dz_plus(k) = (dz(k)+dz(k+1))/2.0
end do
do k = 2, nz
Kvisc_minus(k) = (Kvisc(k)+Kvisc(k-1))/2.0
dz_minus(k) = (dz(k)+dz(k-1))/2.0
end do
do k = 2, nz-1
G(k) = Kvisc_plus(k)/dz_plus(k) * (Field(k+1) - Field(k)) &
- Kvisc_minus(k)/dz_minus(k) * (Field(k) - Field(k-1))
end do
G(1) = Kvisc_plus(1)/dz_plus(1) * (Field(2) - Field(1))
G(nz) = - Kvisc_minus(nz)/dz_minus(nz) * (Field(nz) - Field(nz-1))
do k = 2, nz-1
a(k) = dt / dz(k) * Kvisc_minus(k)/dz_minus(k) * (1.0 - gamma)
b(k) = -1.0 - param%C2_eps * dt * Field(k)/TKE(k) - &
dt / dz(k) * Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma) &
- dt / dz(k) * Kvisc_minus(k)/dz_minus(k) * (1.0 - gamma)
c(k) = dt / dz(k) * Kvisc_plus(k)/dz_plus(k) * (1.0 - gamma)
f(k) = - Field(k) - f_right(k) * dt &
- dt / dz(k) * G(k) * gamma
if (k == 2) then
f(k) = f(k) - a(k) * Field(k-1)
endif
if (k == nz-1) then
f(k) = f(k) - c(k) * Field(k+1)
endif
if (k == 2) a(k) = 0.0
if (k == nz-1) c(k) = 0.0
end do
call tma (Field(2:nz-1), a(2:nz-1), b(2:nz-1), c(2:nz-1), f(2:nz-1), nz - 2)
end subroutine solve_eps_eq_semiimplicit
! --------------------------------------------------------------------------------
subroutine set_config_param(param, tag, ierr)
!< @brief set turbulence closure parameters
! ----------------------------------------------------------------------------
type(kepsilonParamType), intent(out) :: param
integer, intent(out) :: ierr
character(len = *), intent(in) :: tag
integer :: status
! --------------------------------------------------------------------------------
ierr = 0 ! = OK
#ifdef USE_CONFIG_PARSER
call c_config_is_varname(trim(tag)//".C_mu"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".C_mu"//C_NULL_CHAR, param%C_mu, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".C1_eps"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".C1_eps"//C_NULL_CHAR, param%C1_eps, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".C2_eps"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".C2_eps"//C_NULL_CHAR, param%C2_eps, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".C3_eps"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".C3_eps"//C_NULL_CHAR, param%C3_eps, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".ScT_TKE"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".ScT_TKE"//C_NULL_CHAR, param%ScT_TKE, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".ScT_EPS"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".ScT_EPS"//C_NULL_CHAR, param%ScT_EPS, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".TKE_min"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".TKE_min"//C_NULL_CHAR, param%TKE_min, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".EPS_min"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".EPS_min"//C_NULL_CHAR, param%EPS_min, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".kappa"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".kappa"//C_NULL_CHAR, param%kappa, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".PrT0"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".PrT0"//C_NULL_CHAR, param%PrT0, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
call c_config_is_varname(trim(tag)//".c_S2_min"//C_NULL_CHAR, status)
if (status /= 0) then
call c_config_get_float(trim(tag)//".c_S2_min"//C_NULL_CHAR, param%c_S2_min, status)
if (status == 0) then
ierr = 1 ! signal ERROR
return
end if
endif
#else
!> *: just skipping setup without configuration file
#endif
end subroutine
end module