Commit 42345180 authored by Ned Coltman's avatar Ned Coltman
Browse files

[freeflow][rans][sst] Add the SST Turbulence model

parent 9f822132
add_subdirectory(kepsilon)
add_subdirectory(komega)
add_subdirectory(lowrekepsilon)
add_subdirectory(sst)
file(GLOB DUMUX_FREEFLOW_RANS_TWOEQ_HEADERS *.hh *.inc)
install(FILES ${DUMUX_FREEFLOW_RANS_TWOEQ_HEADERS}
......
add_subdirectory(staggered)
file(GLOB DUMUX_FREEFLOW_RANS_TWOEQ_SST_HEADERS *.hh *.inc)
install(FILES ${DUMUX_FREEFLOW_RANS_TWOEQ_SST_HEADERS}
DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/freeflow/rans/twoeq/sst)
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup SSTModel
* \copydoc Dumux::SSTFluxVariables
*/
#ifndef DUMUX_SST_FLUXVARIABLES_HH
#define DUMUX_SST_FLUXVARIABLES_HH
#include <dumux/common/properties.hh>
#include <dumux/freeflow/rans/twoeq/sst/staggered/fluxvariables.hh>
namespace Dumux {
// forward declaration
template<class TypeTag, class BaseFluxVariables, class DiscretizationMethod>
class SSTFluxVariablesImpl;
/*!
* \ingroup SSTModel
* \brief The flux variables class for the SST model.
This is a convenience alias for that actual,
discretization-specific flux variables.
* \note Not all specializations are currently implemented
*/
template<class TypeTag, class BaseFluxVariables>
using SSTFluxVariables = SSTFluxVariablesImpl<TypeTag, BaseFluxVariables, typename GetPropType<TypeTag, Properties::GridGeometry>::DiscretizationMethod>;
} // end namespace
#endif
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup SSTModel
* \copydoc Dumux::SSTIOFields
*/
#ifndef DUMUX_SST_IO_FIELDS_HH
#define DUMUX_SST_IO_FIELDS_HH
#include <dumux/freeflow/rans/iofields.hh>
namespace Dumux {
/*!
* \ingroup SSTModel
* \brief Adds I/O fields for the Reynolds-Averaged Navier-Stokes model
*/
struct SSTIOFields
{
//! Initialize the SSTModel specific output fields.
template <class OutputModule>
static void initOutputModule(OutputModule& out)
{
RANSIOFields::initOutputModule(out);
out.addVolumeVariable([](const auto& v){ return v.turbulentKineticEnergy(); }, "k");
out.addVolumeVariable([](const auto& v){ return v.dissipation(); }, "omega");
}
//! return the names of the primary variables
template <class ModelTraits, class FluidSystem>
static std::string primaryVariableName(int pvIdx = 0, int state = 0)
{
if (pvIdx < ModelTraits::dim() + ModelTraits::numFluidComponents())
return RANSIOFields::template primaryVariableName<ModelTraits, FluidSystem>(pvIdx, state);
else if (pvIdx == ModelTraits::dim() + ModelTraits::numFluidComponents())
return "k";
else
return "omega";
}
};
} // end namespace Dumux
#endif
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup SSTModel
* \copydoc Dumux::SSTResidual
*/
#ifndef DUMUX_SST_LOCAL_RESIDUAL_HH
#define DUMUX_SST_LOCAL_RESIDUAL_HH
#include <dumux/common/properties.hh>
#include <dumux/discretization/method.hh>
#include <dumux/freeflow/navierstokes/localresidual.hh>
#include <dumux/freeflow/rans/twoeq/sst/staggered/localresidual.hh>
namespace Dumux {
// forward declaration
template<class TypeTag, class BaseLocalResidual, class DiscretizationMethod>
class SSTResidualImpl;
/*!
* \ingroup SSTModel
* \brief The local residual class for the SST model.
This is a convenience alias for the actual,
discretization-specific local residual.
* \note Not all specializations are currently implemented
*/
template<class TypeTag, class BaseLocalResidual>
using SSTResidual = SSTResidualImpl<TypeTag, BaseLocalResidual, typename GetPropType<TypeTag, Properties::GridGeometry>::DiscretizationMethod>;
}
#endif
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup SSTModel
*
* \brief A single-phase, isothermal SST 2-Eq. model
*
* \copydoc RANSModel
*
* Two additional PDEs, one for the turbulentKineticEnergy (k) and a second for the dissipation (omega)
* are used to calculate the eddy viscosity for this model.
* The model is taken from Wilcox, 2008 \cite Wilcox2008a.
*
* Turbulent Kinetic Energy balance:
* \f[
* \frac{\partial \varrho k}{\partial t}
* + \nabla \cdot \left( \varrho \boldsymbol{u} k \right)
* - 2\mu_t \boldsymbol{S}\cdot\boldsymbol{S}
* + \beta^* \varrho\omega k
* - \nabla\cdot \left[ \left( \mu + \sigma_k \mu_t \right)\nabla k\right]
* = 0
* \f]
* and \f$ S_{ij} = \frac{1}{2} \left[ \frac{\partial}{\partial x_i} u_j + \frac{\partial}{\partial x_j} u_i \right] \f$
* based on \f$ a_{ij} \cdot b_{ij} = \sum_{i,j} a_{ij} b_{ij} \f$.
*
* Dissipation(rate) balance:
* \f[
* \frac{\partial \varrho \omega}{\partial t}
* + \nabla \cdot \left( \varrho \boldsymbol{u} \omega \right)
* - \frac{\gamma}{\nu_t}\left(2\mu_t \boldsymbol{S}\cdot\boldsymbol{S}\right)
* + \beta^* \varrho\omega^2
* - \nabla\cdot \left[ \left( \mu + \sigma_k \mu_t \right)\nabla \omega\right]
* - 2\varrho\left( 1-F_1\right) \sigma_{\omega 2} \frac{1}{\omega}\nabla k \nabla \omega
* = 0
* \f]
*
* The dynamic eddy viscosity \f$ \mu_\textrm{t} \f$ is calculated as follows:
* \f[ \mu_t = \varrho \frac{a_1 k}{max\left( a_1 \omega; \Omega F_2\right)} \f]
* and \f$ a_1 = 0.31 \f$
* and \f$ \Omega = \sqrt{2\boldsymbol{\Omega}\cdot\boldsymbol{\Omega}} \text{ with } \boldsymbol{\Omega} = \frac{1}{2} \left( \nabla\boldsymbol{u} - \nabla^T\boldsymbol{u}\right) \f$
* and \f$ F_2 = tanh\left( arg_2^2\right) \f$
* and \f$ arg_2 = max\left( 2\frac{\sqrt{k}}{0.09\omega y}; \frac{500\nu}{y^2\omega} \right) \f$
* where y is the distance to the closest wall and \f$ \nu \f$ is the kinematic viscosity.
*/
#ifndef DUMUX_SST_MODEL_HH
#define DUMUX_SST_MODEL_HH
#include <dumux/common/properties.hh>
#include <dumux/freeflow/properties.hh>
#include <dumux/freeflow/rans/model.hh>
#include <dumux/freeflow/rans/twoeq/indices.hh>
#include <dumux/freeflow/turbulencemodel.hh>
#include <dumux/freeflow/rans/problem.hh>
#include "fluxvariables.hh"
#include "localresidual.hh"
#include "volumevariables.hh"
#include "iofields.hh"
namespace Dumux {
namespace Properties {
/*!
*\ingroup SSTModel
* \brief Traits for the sst model
*
* \tparam dimension The dimension of the problem
*/
template<int dimension>
struct SSTModelTraits : RANSModelTraits<dimension>
{
//! The dimension of the model
static constexpr int dim() { return dimension; }
//! There are as many momentum balance equations as dimensions,
//! one mass balance equation and two turbulent transport equations
static constexpr int numEq() { return dim()+1+2; }
//! The number of components
static constexpr int numFluidComponents() { return 1; }
//! The indices
using Indices = RANSTwoEqIndices<dim(), numFluidComponents()>;
//! return the type of turbulence model used
static constexpr auto turbulenceModel()
{ return TurbulenceModel::sst; }
};
///////////////////////////////////////////////////////////////////////////
// default property values for the isothermal SST single phase model
///////////////////////////////////////////////////////////////////////////
// Create new type tags
namespace TTag {
//! The type tag for the single-phase, isothermal SST model
struct SST { using InheritsFrom = std::tuple<RANS>; };
} // end namespace TTag
//! states some specifics of the isothermal SST model
template<class TypeTag>
struct ModelTraits<TypeTag, TTag::SST>
{
private:
using GridView = typename GetPropType<TypeTag, Properties::GridGeometry>::GridView;
static constexpr int dim = GridView::dimension;
public:
using type = SSTModelTraits<dim>;
};
//! The flux variables
template<class TypeTag>
struct FluxVariables<TypeTag, TTag::SST>
{
private:
using BaseFluxVariables = NavierStokesFluxVariables<TypeTag>;
public:
using type = SSTFluxVariables<TypeTag, BaseFluxVariables>;
};
//! The local residual
template<class TypeTag>
struct LocalResidual<TypeTag, TTag::SST>
{
private:
using BaseLocalResidual = NavierStokesResidual<TypeTag>;
public:
using type = SSTResidual<TypeTag, BaseLocalResidual>;
};
//! Set the volume variables property
template<class TypeTag>
struct VolumeVariables<TypeTag, TTag::SST>
{
private:
using PV = GetPropType<TypeTag, Properties::PrimaryVariables>;
using FSY = GetPropType<TypeTag, Properties::FluidSystem>;
using FST = GetPropType<TypeTag, Properties::FluidState>;
using MT = GetPropType<TypeTag, Properties::ModelTraits>;
static_assert(FSY::numPhases == MT::numFluidPhases(), "Number of phases mismatch between model and fluid system");
static_assert(FST::numPhases == MT::numFluidPhases(), "Number of phases mismatch between model and fluid state");
static_assert(!FSY::isMiscible(), "The Navier-Stokes model only works with immiscible fluid systems.");
using Traits = NavierStokesVolumeVariablesTraits<PV, FSY, FST, MT>;
using NSVolVars = NavierStokesVolumeVariables<Traits>;
public:
using type = SSTVolumeVariables<Traits, NSVolVars>;
};
//! The specific I/O fields
template<class TypeTag>
struct IOFields<TypeTag, TTag::SST> { using type = SSTIOFields; };
///////////////////////////////////////////////////////////////////////////
// default property values for the non-isothermal SST single phase model
///////////////////////////////////////////////////////////////////////////
// Create new type tags
namespace TTag {
//! The type tag for the single-phase, non-isothermal SST 2-Eq. model
struct SSTNI { using InheritsFrom = std::tuple<SST, RANSNI>; };
} // end namespace TTag
//! The model traits of the non-isothermal model
template<class TypeTag>
struct ModelTraits<TypeTag, TTag::SSTNI>
{
private:
using GridView = typename GetPropType<TypeTag, Properties::GridGeometry>::GridView;
static constexpr int dim = GridView::dimension;
using IsothermalTraits = SSTModelTraits<dim>;
public:
using type = FreeflowNIModelTraits<IsothermalTraits>;
};
//! Set the volume variables property
template<class TypeTag>
struct VolumeVariables<TypeTag, TTag::SSTNI>
{
private:
using PV = GetPropType<TypeTag, Properties::PrimaryVariables>;
using FSY = GetPropType<TypeTag, Properties::FluidSystem>;
using FST = GetPropType<TypeTag, Properties::FluidState>;
using MT = GetPropType<TypeTag, Properties::ModelTraits>;
static_assert(FSY::numPhases == MT::numFluidPhases(), "Number of phases mismatch between model and fluid system");
static_assert(FST::numPhases == MT::numFluidPhases(), "Number of phases mismatch between model and fluid state");
static_assert(!FSY::isMiscible(), "The Navier-Stokes model only works with immiscible fluid systems.");
using Traits = NavierStokesVolumeVariablesTraits<PV, FSY, FST, MT>;
using NSVolVars = NavierStokesVolumeVariables<Traits>;
public:
using type = SSTVolumeVariables<Traits, NSVolVars>;
};
//! The specific non-isothermal I/O fields
template<class TypeTag>
struct IOFields<TypeTag, TTag::SSTNI> { using type = FreeflowNonIsothermalIOFields<SSTIOFields, true/*turbulenceModel*/>; };
} // end properties
} // end namespace
#endif
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup SSTModel
* \brief SST turbulence model problem base class.
*/
#ifndef DUMUX_SST_PROBLEM_HH
#define DUMUX_SST_PROBLEM_HH
#include <dumux/common/properties.hh>
#include <dumux/common/staggeredfvproblem.hh>
#include <dumux/discretization/localview.hh>
#include <dumux/discretization/staggered/elementsolution.hh>
#include <dumux/discretization/method.hh>
#include <dumux/freeflow/rans/problem.hh>
#include <dumux/freeflow/turbulencemodel.hh>
#include "model.hh"
namespace Dumux {
/*!
* \ingroup SSTModel
* \brief SST turbulence model problem base class.
*
* This implements the 2-equation SST turbulence model developed in Menter1994
*/
template<class TypeTag>
class RANSProblemImpl<TypeTag, TurbulenceModel::sst> : public RANSProblemBase<TypeTag>
{
using ParentType = RANSProblemBase<TypeTag>;
using Implementation = GetPropType<TypeTag, Properties::Problem>;
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
using VolumeVariables = GetPropType<TypeTag, Properties::VolumeVariables>;
using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
using Element = typename GridGeometry::GridView::template Codim<0>::Entity;
using DimVector = typename Element::Geometry::GlobalCoordinate;
public:
RANSProblemImpl(std::shared_ptr<const GridGeometry> gridGeometry, const std::string& paramGroup = "")
: ParentType(gridGeometry, paramGroup)
{ }
/*!
* \brief Correct size of the static (solution independent) wall variables
*/
void updateStaticWallProperties()
{
ParentType::updateStaticWallProperties();
sstModelType_ = getParam<std::string>("RANS.SSTModelType");
// update size and initial values of the global vectors
storedDynamicEddyViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
storedDissipation_.resize(this->gridGeometry().elementMapper().size(), 0.0);
storedDissipationGradient_.resize(this->gridGeometry().elementMapper().size(), DimVector(0.0));
storedTurbulentKineticEnergy_.resize(this->gridGeometry().elementMapper().size(), 0.0);
storedTurbulentKineticEnergyGradient_.resize(this->gridGeometry().elementMapper().size(), DimVector(0.0));
}
/*!
* \brief Update the dynamic (solution dependent) relations to the walls
*
* \param curSol The solution vector.
*/
void updateDynamicWallProperties(const SolutionVector& curSol)
{
ParentType::updateDynamicWallProperties(curSol);
for (const auto& element : elements(this->gridGeometry().gridView()))
{
unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
auto fvGeometry = localView(this->gridGeometry());
fvGeometry.bindElement(element);
for (auto&& scv : scvs(fvGeometry))
{
const int dofIdx = scv.dofIndex();
const auto& cellCenterPriVars = curSol[GridGeometry::cellCenterIdx()][dofIdx];
PrimaryVariables priVars = makePriVarsFromCellCenterPriVars<PrimaryVariables>(cellCenterPriVars);
auto elemSol = elementSolution<typename GridGeometry::LocalView>(std::move(priVars));
// NOTE: first update the turbulence quantities
storedDissipation_[elementIdx] = elemSol[0][Indices::dissipationEqIdx];
storedTurbulentKineticEnergy_[elementIdx] = elemSol[0][Indices::turbulentKineticEnergyEqIdx];
// NOTE: then update the volVars
VolumeVariables volVars;
volVars.update(elemSol, asImp_(), element, scv);
if(sstModelType() == "BSL"){
storedDynamicEddyViscosity_[elementIdx] = volVars.calculateEddyViscosityBSL();
}
else if(sstModelType() == "SST"){
storedDynamicEddyViscosity_[elementIdx] = volVars.calculateEddyViscositySST();
}
}
}
// calculate cell-centered gradients
for (const auto& element : elements(this->gridGeometry().gridView()))
{
const unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
for (unsigned int dimIdx = 0; dimIdx < DimVector::dimension; ++dimIdx)
{
const unsigned neighborIdx0 = ParentType::neighborIndex(elementIdx, dimIdx, 0);
const unsigned neighborIdx1 = ParentType::neighborIndex(elementIdx, dimIdx, 1);
// Cell centered TKE Gradient
storedTurbulentKineticEnergyGradient_[elementIdx][dimIdx]
= (storedTurbulentKineticEnergy(neighborIdx1) - storedTurbulentKineticEnergy(neighborIdx0))
/ (ParentType::cellCenter(neighborIdx1)[dimIdx] - ParentType::cellCenter(neighborIdx0)[dimIdx]);
// Cell centered Omega Gradient
storedDissipationGradient_[elementIdx][dimIdx]
= (storedDissipation(neighborIdx1) - storedDissipation(neighborIdx0))
/ (ParentType::cellCenter(neighborIdx1)[dimIdx] - ParentType::cellCenter(neighborIdx0)[dimIdx]);
}
}
}
//! \brief Returns the \f$ \beta_{\omega} \f$ constant
const Scalar betaOmega() const
{ return 0.0708; }
bool useStoredEddyViscosity() const
{
static const bool useStoredEddyViscosity = getParamFromGroup<bool>(this->paramGroup(), "RANS.UseStoredEddyViscosity", false);
return useStoredEddyViscosity;
}
Scalar storedDynamicEddyViscosity(const int elementIdx) const
{ return storedDynamicEddyViscosity_[elementIdx]; }
Scalar storedTurbulentKineticEnergy(const int elementIdx) const
{ return storedTurbulentKineticEnergy_[elementIdx]; }
Scalar storedDissipation(const int elementIdx) const
{ return storedDissipation_[elementIdx]; }
DimVector storedTurbulentKineticEnergyGradient(const int elementIdx) const
{ return storedTurbulentKineticEnergyGradient_[elementIdx]; }
DimVector storedDissipationGradient(const int elementIdx) const
{ return storedDissipationGradient_[elementIdx]; }
std::string sstModelType() const
{ return sstModelType_; }
private:
std::vector<Scalar> storedDynamicEddyViscosity_;
std::vector<Scalar> storedTurbulentKineticEnergy_;
std::vector<Scalar> storedDissipation_;
std::vector<DimVector> storedDissipationGradient_;
std::vector<DimVector> storedTurbulentKineticEnergyGradient_;
std::string sstModelType_;
//! Returns the implementation of the problem (i.e. static polymorphism)
Implementation &asImp_()
{ return *static_cast<Implementation *>(this); }
//! \copydoc asImp_()
const Implementation &asImp_() const
{ return *static_cast<const Implementation *>(this); }
};
} // end namespace Dumux
#endif
file(GLOB DUMUX_FREEFLOW_RANS_TWOEQ_SST_STAGGERED_HEADERS *.hh *.inc)
install(FILES ${DUMUX_FREEFLOW_RANS_TWOEQ_SST_STAGGERED_HEADERS}
DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/freeflow/rans/twoeq/sst/staggered)