diff --git a/obl_k_epsilon.f90 b/obl_k_epsilon.f90
index 182dd0642a652381936bb84a36e1bfccfc924017..f81d0259ce21dba2b1418092053e311eabd548cc 100644
--- a/obl_k_epsilon.f90
+++ b/obl_k_epsilon.f90
@@ -4,6 +4,11 @@ module obl_k_epsilon
! 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
@@ -25,20 +30,24 @@ module obl_k_epsilon
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 :: PrT0 = 1.25 !< turbulent Prandtl number
- 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 :: 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 :: Sch_TKE = 1.0 !< Schmidt number for TKE
- real :: Sch_eps = 1.11 !< Scmidt number for eps
+ 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
@@ -270,7 +279,7 @@ module obl_k_epsilon
real :: B_TKE(nz) !< buoyancy production of TKE, [m**2 / s**3]
real :: eps(nz) !< TKE dissipation rate, [W/kg]
- Kvisc(:) = Km(:) / param%Sch_TKE
+ Kvisc(:) = Km(:) / param%ScT_TKE
gamma = 0.0
@@ -344,7 +353,7 @@ module obl_k_epsilon
real :: eps(nz) !< TKE dissipation rate, [W/kg]
real :: TKE(nz) !< TKE, [J/kg]
- Kvisc(:) = Km(:) / param%Sch_eps
+ Kvisc(:) = Km(:) / param%ScT_EPS
gamma = 0.0
@@ -416,7 +425,7 @@ module obl_k_epsilon
real :: B_TKE(nz) !< buoyancy production of TKE, [m**2 / s**3]
real :: eps(nz) !< TKE dissipation rate, [W/kg]
- Kvisc(:) = Km(:) / param%Sch_TKE
+ Kvisc(:) = Km(:) / param%ScT_TKE
gamma = 0.0
@@ -491,7 +500,7 @@ module obl_k_epsilon
real :: TKE(nz) !< TKE, [J/kg]
- Kvisc(:) = Km(:) / param%Sch_eps
+ Kvisc(:) = Km(:) / param%ScT_EPS
gamma = 0.0
@@ -540,4 +549,113 @@ module obl_k_epsilon
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
\ No newline at end of file
diff --git a/obl_main.f90 b/obl_main.f90
index 8c7fcc9fe60453708b2d6d081cd8e31a3ce0b53f..dbc36936264e2d23e5ee19fe7e9c95aed2cb447e 100644
--- a/obl_main.f90
+++ b/obl_main.f90
@@ -26,14 +26,17 @@ program obl_main
use obl_diag
use obl_bc
use obl_pph, only: pphParamType, &
+ set_config_param_pph => set_config_param, &
define_stability_functions_pph => define_stability_functions, &
init_turbulence_closure_pph => init_turbulence_closure, &
advance_turbulence_closure_pph => advance_turbulence_closure
use obl_pph_dyn, only: pphDynParamType, &
+ set_config_param_pph_dyn => set_config_param, &
define_stability_functions_pph_dyn => define_stability_functions, &
init_turbulence_closure_pph_dyn => init_turbulence_closure, &
advance_turbulence_closure_pph_dyn => advance_turbulence_closure
use obl_k_epsilon, only: kepsilonParamType, &
+ set_config_param_k_epsilon => set_config_param, &
define_stability_functions_k_epsilon => define_stability_functions, &
init_turbulence_closure_k_epsilon => init_turbulence_closure, &
advance_turbulence_closure_k_epsilon => advance_turbulence_closure, &
@@ -238,6 +241,25 @@ program obl_main
endif
! ----------------------------------------------------------------------------
+ !< setting turbulence closure parameters
+ ! ----------------------------------------------------------------------------
+ call set_config_param_pph(param_pph, "turbulence_model.pph", ierr)
+ if (ierr /= 0) then
+ write(*, *) ' FAILURE! > unable to set ''pph'' parameters '
+ return
+ endif
+ call set_config_param_pph_dyn(param_pph_dyn, "turbulence_model.pph_dyn", ierr)
+ if (ierr /= 0) then
+ write(*, *) ' FAILURE! > unable to set ''pph-dyn'' parameters '
+ return
+ endif
+ call set_config_param_k_epsilon(param_k_epsilon, "turbulence_model.k_epsilon", ierr)
+ if (ierr /= 0) then
+ write(*, *) ' FAILURE! > unable to set ''k-epsilon'' parameters '
+ return
+ endif
+ ! ----------------------------------------------------------------------------
+
!< initialization of main fields
! ----------------------------------------------------------------------------
call set_initial_conditions(grid, obl_config_id, ierr)
diff --git a/obl_pph.f90 b/obl_pph.f90
index fa63efb63f3d49b6915b4898bf758c6cb612fbb9..74c3d75ef97bca832bdcb5a10d1c98b4acbe5702 100644
--- a/obl_pph.f90
+++ b/obl_pph.f90
@@ -1,9 +1,16 @@
+#include "obl_def.fi"
+
module obl_pph
!< @brief standard pacanowski-philander scheme
! --------------------------------------------------------------------------------
! modules used
! --------------------------------------------------------------------------------
+#ifdef USE_CONFIG_PARSER
+ use iso_c_binding, only: C_NULL_CHAR
+ use config_parser
+#endif
+
use obl_grid
use obl_bc
use obl_turb_common
@@ -21,6 +28,8 @@ module obl_pph
public :: get_eddy_viscosity, get_eddy_diffusivity
public :: get_eddy_viscosity_vec, get_eddy_diffusivity_vec
+
+ public :: set_config_param
! --------------------------------------------------------------------------------
!> @brief Pacanowski-Philander parameters
@@ -154,5 +163,98 @@ module obl_pph
Kh(k) = param%Kh_0 / (1.0 + param%alpha * Ri_grad(k))**(param%n + 1.0) + param%Kh_b
end do
end subroutine
+
+ ! --------------------------------------------------------------------------------
+ subroutine set_config_param(param, tag, ierr)
+ !< @brief set turbulence closure parameters
+ ! ----------------------------------------------------------------------------
+ type(pphParamType), 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)//".Km_0"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".Km_0"//C_NULL_CHAR, param%Km_0, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".Kh_0"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".Kh_0"//C_NULL_CHAR, param%Kh_0, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".alpha"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".alpha"//C_NULL_CHAR, param%alpha, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".n"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".n"//C_NULL_CHAR, param%n, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".Km_b"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".Km_b"//C_NULL_CHAR, param%Km_b, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".Kh_b"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".Kh_b"//C_NULL_CHAR, param%Kh_b, 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
diff --git a/obl_pph_dyn.f90 b/obl_pph_dyn.f90
index ebb94268382557a4a4a54a84c35196b3662f1458..2c141a4794827529da87d1a421112da3acf6ba61 100644
--- a/obl_pph_dyn.f90
+++ b/obl_pph_dyn.f90
@@ -5,6 +5,11 @@ module obl_pph_dyn
! modules used
! --------------------------------------------------------------------------------
+#ifdef USE_CONFIG_PARSER
+ use iso_c_binding, only: C_NULL_CHAR
+ use config_parser
+#endif
+
use obl_grid
use obl_bc
use obl_turb_common
@@ -23,6 +28,8 @@ module obl_pph_dyn
public :: get_eddy_viscosity, get_eddy_diffusivity
public :: get_eddy_viscosity_vec, get_eddy_diffusivity_vec
+
+ public :: set_config_param
! --------------------------------------------------------------------------------
!> @brief Pacanowski-Philander dynamic parameters
@@ -32,11 +39,11 @@ module obl_pph_dyn
real :: Km_b = 0.0001 !< background eddy viscosity, [m**2 / s], *: INMCM: 5 * 10**(-6)
real :: Kh_b = 0.00001 !< background eddy diffusivity, [m**2 / s]
+ real :: C_l = 4.0 / 9.0 !< length scale tuning constant, [n/d], *: = (4.0 / pi^2) in Obukhov length scale
+
real :: kappa = 0.4 !< von Karman constant
real :: PrT0 = 1.0 !< turbulent Prandtl number, *: probably redundant
- real :: C_l = 4.0 / 9.0 !< length scale tuning constant, [n/d], *: = (4.0 / pi^2) in Obukhov length scale
-
real :: c_S2_min = 1e-5 !< min shear frequency, [(1/s)**2]
end type
@@ -187,4 +194,89 @@ module obl_pph_dyn
end do
end subroutine
+ ! --------------------------------------------------------------------------------
+ subroutine set_config_param(param, tag, ierr)
+ !< @brief set turbulence closure parameters
+ ! ----------------------------------------------------------------------------
+ type(pphDynParamType), 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)//".alpha"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".alpha"//C_NULL_CHAR, param%alpha, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".n"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".n"//C_NULL_CHAR, param%n, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".Km_b"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".Km_b"//C_NULL_CHAR, param%Km_b, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".Kh_b"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".Kh_b"//C_NULL_CHAR, param%Kh_b, status)
+ if (status == 0) then
+ ierr = 1 ! signal ERROR
+ return
+ end if
+ endif
+ call c_config_is_varname(trim(tag)//".C_l"//C_NULL_CHAR, status)
+ if (status /= 0) then
+ call c_config_get_float(trim(tag)//".C_l"//C_NULL_CHAR, param%C_l, 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
\ No newline at end of file