diff --git a/test/freeflow/rans/problem.hh b/test/freeflow/rans/problem.hh index 4bfef43a545993089f5807d0214844ce9bf75a95..6245a4d044e72f2b40b0a035d83383e35613b406 100644 --- a/test/freeflow/rans/problem.hh +++ b/test/freeflow/rans/problem.hh @@ -272,31 +272,17 @@ public: * \param scv the sub control volume * \note used for cell-centered discretization schemes */ - template<bool enable = (ModelTraits::turbulenceModel() == TurbulenceModel::komega - || ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon), - std::enable_if_t<!enable, int> = 0> PrimaryVariables dirichlet(const Element& element, const SubControlVolume& scv) const { - const auto globalPos = scv.center(); - PrimaryVariables values(initialAtPos(globalPos)); - return values; - } - - /*! - * \brief Evaluate the boundary conditions for fixed values at cell centers - * - * \param element The finite element - * \param scv the sub control volume - * \note used for cell-centered discretization schemes - */ - template<bool enable = (ModelTraits::turbulenceModel() == TurbulenceModel::komega - || ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon), - std::enable_if_t<enable, int> = 0> - PrimaryVariables dirichlet(const Element& element, const SubControlVolume& scv) const - { - using SetDirichletCellForBothTurbEq = std::integral_constant<bool, (ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon)>; - - return dirichletTurbulentTwoEq_(element, scv, SetDirichletCellForBothTurbEq{}); + if constexpr (ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon + || ModelTraits::turbulenceModel() == TurbulenceModel::komega) + return dirichletTurbulentTwoEq_(element, scv); + else + { + const auto globalPos = scv.center(); + PrimaryVariables values(initialAtPos(globalPos)); + return values; + } } /*! @@ -428,21 +414,32 @@ private: return ParentType::isDirichletCell(element, fvGeometry, scv, pvIdx); } - //! Specialization for the KOmega + //! Specialization for the kepsilon and komega template<class Element, class SubControlVolume> PrimaryVariables dirichletTurbulentTwoEq_(const Element& element, - const SubControlVolume& scv, - std::false_type) const + const SubControlVolume& scv) const { const auto globalPos = scv.center(); PrimaryVariables values(initialAtPos(globalPos)); unsigned int elementIdx = this->gridGeometry().elementMapper().index(element); - const auto wallDistance = ParentType::wallDistance_[elementIdx]; - using std::pow; - values[Indices::dissipationEqIdx] = 6.0 * ParentType::kinematicViscosity_[elementIdx] - / (ParentType::betaOmega() * pow(wallDistance, 2)); - return values; + if constexpr (ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon) + { + // For the kepsilon model we set a fixed value for the turbulent kinetic energy and the dissipation + values[Indices::turbulentKineticEnergyEqIdx] = this->turbulentKineticEnergyWallFunction(elementIdx); + values[Indices::dissipationEqIdx] = this->dissipationWallFunction(elementIdx); + return values; + } + else + { + static_assert(ModelTraits::turbulenceModel() == TurbulenceModel::komega, "Only valid for Komega"); + // For the komega model we set a fixed value for the dissipation + const auto wallDistance = ParentType::wallDistance_[elementIdx]; + using std::pow; + values[Indices::dissipationEqIdx] = 6.0 * ParentType::kinematicViscosity_[elementIdx] + / (ParentType::betaOmega() * wallDistance * wallDistance); + return values; + } } Scalar eps_; diff --git a/test/freeflow/ransnc/problem.hh b/test/freeflow/ransnc/problem.hh index 530776ea87acc4efece964a31d0c59bd5922eaa4..92d41877b93207defdd7caf5243a7be013524031 100644 --- a/test/freeflow/ransnc/problem.hh +++ b/test/freeflow/ransnc/problem.hh @@ -253,11 +253,7 @@ public: const FVElementGeometry& fvGeometry, const SubControlVolume& scv, int pvIdx) const - { - using IsKOmegaKEpsilon = std::integral_constant<bool, (ModelTraits::turbulenceModel() == TurbulenceModel::komega - || ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon)>; - return isDirichletCell_(element, fvGeometry, scv, pvIdx, IsKOmegaKEpsilon{}); - } + { return isDirichletCell_(element, fvGeometry, scv, pvIdx); } /*! * \brief Evaluate the boundary conditions for a dirichlet values at the boundary. @@ -295,31 +291,17 @@ public: * \param scv the sub control volume * \note used for cell-centered discretization schemes */ - template<bool enable = (ModelTraits::turbulenceModel() == TurbulenceModel::komega - || ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon), - std::enable_if_t<!enable, int> = 0> - PrimaryVariables dirichlet(const Element& element, const SubControlVolume& scv) const - { - const auto globalPos = scv.center(); - PrimaryVariables values(initialAtPos(globalPos)); - return values; - } - - /*! - * \brief Evaluate the boundary conditions for fixed values at cell centers - * - * \param element The finite element - * \param scv the sub control volume - * \note used for cell-centered discretization schemes - */ - template<bool enable = (ModelTraits::turbulenceModel() == TurbulenceModel::komega - || ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon), - std::enable_if_t<enable, int> = 0> PrimaryVariables dirichlet(const Element& element, const SubControlVolume& scv) const { - using SetDirichletCellForBothTurbEq = std::integral_constant<bool, (ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon)>; - - return dirichletTurbulentTwoEq_(element, scv, SetDirichletCellForBothTurbEq{}); + if constexpr (ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon + || ModelTraits::turbulenceModel() == TurbulenceModel::komega) + return dirichletTurbulentTwoEq_(element, scv); + else + { + const auto globalPos = scv.center(); + PrimaryVariables values(initialAtPos(globalPos)); + return values; + } } /*! @@ -423,67 +405,32 @@ private: } } - //! Forward to ParentType template<class Element, class FVElementGeometry, class SubControlVolume> bool isDirichletCell_(const Element& element, const FVElementGeometry& fvGeometry, const SubControlVolume& scv, - int pvIdx, - std::false_type) const - { - return ParentType::isDirichletCell(element, fvGeometry, scv, pvIdx); - } - - //! Specialization for the KOmega and KEpsilon Models - template<class Element, class FVElementGeometry, class SubControlVolume> - bool isDirichletCell_(const Element& element, - const FVElementGeometry& fvGeometry, - const SubControlVolume& scv, - int pvIdx, - std::true_type) const - { - using SetDirichletCellForBothTurbEq = std::integral_constant<bool, (ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon)>; - return isDirichletCellTurbulentTwoEq_(element, fvGeometry, scv, pvIdx, SetDirichletCellForBothTurbEq{}); - } - - //! Specialization for the KEpsilon Model - template<class Element> - bool isDirichletCellTurbulentTwoEq_(const Element& element, - const FVElementGeometry& fvGeometry, - const SubControlVolume& scv, - int pvIdx, - std::true_type) const - { - const auto eIdx = this->gridGeometry().elementMapper().index(element); - - // set a fixed turbulent kinetic energy and dissipation near the wall - if (this->inNearWallRegion(eIdx)) - return pvIdx == Indices::turbulentKineticEnergyEqIdx || pvIdx == Indices::dissipationEqIdx; - - // set a fixed dissipation at the matching point - if (this->isMatchingPoint(eIdx)) - return pvIdx == Indices::dissipationEqIdx;// set a fixed dissipation (omega) for all cells at the wall - - return false; - } - - //! Specialization for the KOmega Model - template<class Element> - bool isDirichletCellTurbulentTwoEq_(const Element& element, - const FVElementGeometry& fvGeometry, - const SubControlVolume& scv, - int pvIdx, - std::false_type) const - { - // set a fixed dissipation (omega) for all cells at the wall - for (const auto& scvf : scvfs(fvGeometry)) - if (isOnWallAtPos(scvf.center()) && pvIdx == Indices::dissipationIdx) - return true; - - return false; - - } + const int& pvIdx) const { + if constexpr (ModelTraits::turbulenceModel() == TurbulenceModel::kepsilon) + { + const auto eIdx = this->gridGeometry().elementMapper().index(element); + // For the kepsilon model we set fixed values within the matching point and at the wall + if (this->inNearWallRegion(eIdx)) + return pvIdx == Indices::turbulentKineticEnergyEqIdx || pvIdx == Indices::dissipationEqIdx; + if (this->isMatchingPoint(eIdx)) + return pvIdx == Indices::dissipationEqIdx; + return false; + } + else if constexpr (ModelTraits::turbulenceModel() == TurbulenceModel::komega) + { + // For the komega model we set a fixed dissipation (omega) for all cells at the wall + for (const auto& scvf : scvfs(fvGeometry)) + if (isOnWallAtPos(scvf.center()) && pvIdx == Indices::dissipationIdx) + return true; + return false; + } + else + return ParentType::isDirichletCell(element, fvGeometry, scv, pvIdx); } //! Specialization for the kepsilon and komega @@ -509,7 +456,7 @@ private: const auto wallDistance = ParentType::wallDistance_[elementIdx]; using std::pow; values[Indices::dissipationEqIdx] = 6.0 * ParentType::kinematicViscosity_[elementIdx] - / (ParentType::betaOmega() * pow(wallDistance, 2)); + / (ParentType::betaOmega() * wallDistance * wallDistance); return values; } }