25#ifndef DUMUX_STOKES1P2C_SUBPROBLEM_HH
26#define DUMUX_STOKES1P2C_SUBPROBLEM_HH
28#include <dune/grid/yaspgrid.hh>
39template <
class TypeTag>
40class StokesSubProblem;
54template<
class TypeTag>
55struct Grid<TypeTag, TTag::StokesOnePTwoC> {
using type = Dune::YaspGrid<2, Dune::TensorProductCoordinates<GetPropType<TypeTag, Properties::Scalar>, 2> >; };
58template<
class TypeTag>
62 static constexpr auto phaseIdx = H2OAir::gasPhaseIdx;
66template<
class TypeTag>
67struct ReplaceCompEqIdx<TypeTag, TTag::StokesOnePTwoC> {
static constexpr int value = 3; };
70template<
class TypeTag>
71struct UseMoles<TypeTag, TTag::StokesOnePTwoC> {
static constexpr bool value =
true; };
74template<
class TypeTag>
77template<
class TypeTag>
79template<
class TypeTag>
80struct EnableGridFluxVariablesCache<TypeTag, TTag::StokesOnePTwoC> {
static constexpr bool value =
true; };
81template<
class TypeTag>
82struct EnableGridVolumeVariablesCache<TypeTag, TTag::StokesOnePTwoC> {
static constexpr bool value =
true; };
91template <
class TypeTag>
92class StokesSubProblem :
public NavierStokesProblem<TypeTag>
94 using ParentType = NavierStokesProblem<TypeTag>;
96 using GridView = GetPropType<TypeTag, Properties::GridView>;
97 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
98 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
100 using BoundaryTypes = GetPropType<TypeTag, Properties::BoundaryTypes>;
102 using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
103 using FVElementGeometry =
typename GridGeometry::LocalView;
104 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
105 using Element =
typename GridView::template Codim<0>::Entity;
106 using ElementVolumeVariables =
typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView;
107 using ElementFaceVariables =
typename GetPropType<TypeTag, Properties::GridFaceVariables>::LocalView;
108 using FluidState = GetPropType<TypeTag, Properties::FluidState>;
110 using GlobalPosition =
typename Element::Geometry::GlobalCoordinate;
112 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
113 using NumEqVector = GetPropType<TypeTag, Properties::NumEqVector>;
115 using CouplingManager = GetPropType<TypeTag, Properties::CouplingManager>;
116 using TimeLoopPtr = std::shared_ptr<TimeLoop<Scalar>>;
120 static constexpr bool useMoles = GetPropType<TypeTag, Properties::ModelTraits>::useMoles();
124 : ParentType(gridGeometry,
"Stokes"), eps_(1e-6), couplingManager_(
couplingManager)
126 refVelocity_ = getParamFromGroup<Scalar>(this->paramGroup(),
"Problem.RefVelocity");
127 refPressure_ = getParamFromGroup<Scalar>(this->paramGroup(),
"Problem.RefPressure");
128 refMoleFrac_ = getParamFromGroup<Scalar>(this->paramGroup(),
"Problem.refMoleFrac");
129 refTemperature_ = getParamFromGroup<Scalar>(this->paramGroup(),
"Problem.RefTemperature");
132 getParamFromGroup<std::string>(this->paramGroup(),
"Problem.InterfaceDiffusionCoefficientAvg"));
133 problemName_ = getParam<std::string>(
"Vtk.OutputName") +
"_" + getParamFromGroup<std::string>(this->paramGroup(),
"Problem.Name");
139 const std::string&
name()
const
153 {
return refTemperature_; }
161 {
return NumEqVector(0.0); }
177 const SubControlVolumeFace& scvf)
const
179 BoundaryTypes values;
181 const auto& globalPos = scvf.center();
184 values.setNeumann(Indices::energyEqIdx);
187 if (onUpperBoundary_(globalPos) || onLeftBoundary_(globalPos))
189 values.setDirichlet(Indices::velocityXIdx);
190 values.setDirichlet(Indices::velocityYIdx);
191 values.setNeumann(Indices::conti0EqIdx);
192 values.setNeumann(Indices::conti0EqIdx + 1);
195 if (onRightBoundary_(globalPos))
197 values.setDirichlet(Indices::pressureIdx);
198 values.setOutflow(Indices::conti0EqIdx + 1);
201 values.setOutflow(Indices::energyEqIdx);
205 if (
couplingManager().isCoupledEntity(CouplingManager::stokesIdx, scvf))
207 values.setCouplingNeumann(Indices::conti0EqIdx);
208 values.setCouplingNeumann(Indices::conti0EqIdx + 1);
209 values.setCouplingNeumann(Indices::momentumYBalanceIdx);
210 values.setBJS(Indices::momentumXBalanceIdx);
220 PrimaryVariables values(0.0);
236 const FVElementGeometry& fvGeometry,
237 const ElementVolumeVariables& elemVolVars,
238 const ElementFaceVariables& elemFaceVars,
239 const SubControlVolumeFace& scvf)
const
241 PrimaryVariables values(0.0);
242 const auto& globalPos = scvf.dofPosition();
243 const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
245 FluidState fluidState;
246 updateFluidStateForBC_(fluidState, elemVolVars[scv].
pressure());
248 const Scalar
density = useMoles ? fluidState.molarDensity(0) : fluidState.density(0);
249 const Scalar xVelocity = xVelocity_(globalPos);
251 if (onLeftBoundary_(globalPos))
258 values[Indices::energyEqIdx] = -xVelocity * fluidState.density(0) * fluidState.enthalpy(0);
262 if(
couplingManager().isCoupledEntity(CouplingManager::stokesIdx, scvf))
264 values[Indices::momentumYBalanceIdx] =
couplingManager().couplingData().momentumCouplingCondition(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
266 const auto massFlux =
couplingManager().couplingData().massCouplingCondition(element, fvGeometry, elemVolVars, elemFaceVars, scvf, diffCoeffAvgType_);
267 values[Indices::conti0EqIdx] = massFlux[0];
268 values[Indices::conti0EqIdx + 1] = massFlux[1];
271 values[Indices::energyEqIdx] =
couplingManager().couplingData().energyCouplingCondition(element, fvGeometry, elemVolVars, elemFaceVars, scvf, diffCoeffAvgType_);
282 {
return *couplingManager_; }
296 FluidState fluidState;
301 PrimaryVariables values(0.0);
302 values[Indices::pressureIdx] =
refPressure() +
density*this->
gravity()[1]*(globalPos[1] - this->gridGeometry().bBoxMin()[1]);
304 values[Indices::velocityXIdx] = xVelocity_(globalPos);
315 {
return refVelocity_ ;}
319 {
return refPressure_; }
323 {
return refMoleFrac_; }
327 {
return refTemperature_; }
331 { timeLoop_ = timeLoop; }
334 {
return timeLoop_->time(); }
340 Scalar
permeability(
const Element& element,
const SubControlVolumeFace& scvf)
const
342 return couplingManager().couplingData().darcyPermeability(element, scvf);
349 Scalar
alphaBJ(
const SubControlVolumeFace& scvf)
const
351 return couplingManager().problem(CouplingManager::darcyIdx).spatialParams().beaversJosephCoeffAtPos(scvf.center());
357 bool onLeftBoundary_(
const GlobalPosition &globalPos)
const
358 {
return globalPos[0] < this->gridGeometry().bBoxMin()[0] + eps_; }
360 bool onRightBoundary_(
const GlobalPosition &globalPos)
const
361 {
return globalPos[0] > this->gridGeometry().bBoxMax()[0] - eps_; }
363 bool onLowerBoundary_(
const GlobalPosition &globalPos)
const
364 {
return globalPos[1] < this->gridGeometry().bBoxMin()[1] + eps_; }
366 bool onUpperBoundary_(
const GlobalPosition &globalPos)
const
367 {
return globalPos[1] > this->gridGeometry().bBoxMax()[1] - eps_; }
370 void updateFluidStateForBC_(FluidState& fluidState,
const Scalar
pressure)
const
373 fluidState.setPressure(0,
pressure);
374 fluidState.setSaturation(0, 1.0);
376 fluidState.setMoleFraction(0, 0, 1.0 -
refMoleFrac());
378 typename FluidSystem::ParameterCache paramCache;
379 paramCache.updatePhase(fluidState, 0);
382 fluidState.setDensity(0,
density);
387 const Scalar enthalpy = FluidSystem::enthalpy(fluidState, paramCache, 0);
388 fluidState.setEnthalpy(0, enthalpy);
392 const Scalar xVelocity_(
const GlobalPosition &globalPos)
const
395 return 4 * vmax * (globalPos[1] - this->gridGeometry().bBoxMin()[1]) * (this->gridGeometry().bBoxMax()[1] - globalPos[1])
396 / (height_() * height_());
400 const Scalar height_()
const
401 {
return this->gridGeometry().bBoxMax()[1] - this->gridGeometry().bBoxMin()[1]; }
408 Scalar refTemperature_;
409 std::string problemName_;
410 TimeLoopPtr timeLoop_;
412 std::shared_ptr<CouplingManager> couplingManager_;
414 DiffusionCoefficientAveragingType diffCoeffAvgType_;
An adapter for multi-phase fluid systems to be used with (compositional) one-phase models.
A compositional two-phase fluid system with water and air as components in both, the liquid and the g...
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
Property tag Indices
Definition: porousmediumflow/sequential/properties.hh:59
std::string molarDensity(int phaseIdx) noexcept
I/O name of molar density for multiphase systems.
Definition: name.hh:83
std::string pressure(int phaseIdx) noexcept
I/O name of pressure for multiphase systems.
Definition: name.hh:34
std::string density(int phaseIdx) noexcept
I/O name of density for multiphase systems.
Definition: name.hh:65
The DUNE grid type.
Definition: common/properties.hh:57
Property to specify the type of a problem which has to be solved.
Definition: common/properties.hh:69
The component balance index that should be replaced by the total mass/mole balance.
Definition: common/properties.hh:104
Definition: common/properties.hh:169
The type of the fluid system to use.
Definition: common/properties.hh:223
const GravityVector & gravity() const
Returns the acceleration due to gravity.
Definition: dumux/freeflow/navierstokes/problem.hh:134
An adapter for multi-phase fluid systems to be used with (compositional) one-phase models.
Definition: 1padapter.hh:46
A compositional two-phase fluid system with water and air as components in both, the liquid and the g...
Definition: h2oair.hh:75
DiffusionCoefficientAveragingType
Defines which kind of averanging of diffusion coefficiencients (moleculat diffusion or thermal conduc...
Definition: couplingdata.hh:52
static DiffusionCoefficientAveragingType stringToEnum(DiffusionCoefficientAveragingType, const std::string &diffusionCoefficientAveragingType)
Convenience function to convert user input given as std::string to the corresponding enum class used ...
Definition: couplingdata.hh:60
Test problem for the 1pnc (Navier-) Stokes problem.
Definition: 1p_2p/problem_stokes.hh:68
Scalar temperature() const
Returns the temperature within the domain in [K].
Definition: 1p2c_2p2c/problem_stokes.hh:152
const Scalar refPressure() const
Returns the reference pressure.
Definition: 1p2c_2p2c/problem_stokes.hh:318
PrimaryVariables dirichletAtPos(const GlobalPosition &pos) const
Evaluates the boundary conditions for a Dirichlet control volume.
Definition: 1p2c_2p2c/problem_stokes.hh:218
void setTimeLoop(TimeLoopPtr timeLoop)
Definition: 1p2c_2p2c/problem_stokes.hh:330
Scalar alphaBJ(const SubControlVolumeFace &scvf) const
Returns the alpha value required as input parameter for the Beavers-Joseph-Saffman boundary condition...
Definition: 1p2c_2p2c/problem_stokes.hh:349
const Scalar refMoleFrac() const
Returns the reference mass fraction.
Definition: 1p2c_2p2c/problem_stokes.hh:322
BoundaryTypes boundaryTypes(const Element &element, const SubControlVolumeFace &scvf) const
Specifies which kind of boundary condition should be used for which equation on a given boundary segm...
Definition: 1p2c_2p2c/problem_stokes.hh:176
const CouplingManager & couplingManager() const
Get the coupling manager.
Definition: 1p2c_1p2c/diffusionlawcomparison/problem_stokes.hh:254
StokesSubProblem(std::shared_ptr< const GridGeometry > gridGeometry, std::shared_ptr< CouplingManager > couplingManager)
Definition: 1p2c_2p2c/problem_stokes.hh:123
NumEqVector sourceAtPos(const GlobalPosition &globalPos) const
Returns the sources within the domain.
Definition: 1p2c_2p2c/problem_stokes.hh:160
NumEqVector neumann(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf) const
Evaluates the boundary conditions for a Neumann control volume.
Definition: 1p2c_2p2c/problem_stokes.hh:235
const Scalar refTemperature() const
Returns the reference temperature.
Definition: 1p2c_2p2c/problem_stokes.hh:326
Scalar time() const
Definition: 1p2c_2p2c/problem_stokes.hh:333
const Scalar refVelocity() const
Returns the reference velocity.
Definition: 1p2c_2p2c/problem_stokes.hh:314
const std::string & name() const
The problem name.
Definition: 1p2c_2p2c/problem_stokes.hh:139
Scalar permeability(const Element &element, const SubControlVolumeFace &scvf) const
Returns the intrinsic permeability of required as input parameter for the Beavers-Joseph-Saffman boun...
Definition: 1p2c_2p2c/problem_stokes.hh:340
PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
Evaluates the initial value for a control volume.
Definition: 1p2c_1p2c/diffusionlawcomparison/problem_stokes.hh:267
Definition: 1p2c_1p2c/diffusionlawcomparison/problem_stokes.hh:57
std::tuple< NavierStokesNC, StaggeredFreeFlowModel > InheritsFrom
Definition: 1p2c_1p2c/diffusionlawcomparison/problem_stokes.hh:57
Dune::YaspGrid< 2, Dune::EquidistantOffsetCoordinates< GetPropType< TypeTag, Properties::Scalar >, 2 > > type
Definition: 1p2c_1p2c/diffusionlawcomparison/problem_stokes.hh:71
Defines a type tag and some properties for ree-flow models using the staggered scheme....