25#ifndef DUMUX_CHANNEL_NC_TEST_PROBLEM_HH
26#define DUMUX_CHANNEL_NC_TEST_PROBLEM_HH
29#define ENABLECACHING 1
32#include <dune/grid/yaspgrid.hh>
45template <
class TypeTag>
46class ChannelNCTestProblem;
60template<
class TypeTag>
64 static constexpr int phaseIdx = H2OAir::liquidPhaseIdx;
68template<
class TypeTag>
69struct ReplaceCompEqIdx<TypeTag, TTag::ChannelNCTest> {
static constexpr int value = 0; };
72template<
class TypeTag>
73struct Grid<TypeTag, TTag::ChannelNCTest> {
using type = Dune::YaspGrid<2>; };
76template<
class TypeTag>
79template<
class TypeTag>
81template<
class TypeTag>
83template<
class TypeTag>
85template<
class TypeTag>
90template<
class TypeTag>
91struct UseMoles<TypeTag, TTag::ChannelNCTest> {
static constexpr bool value =
false; };
93template<
class TypeTag>
94struct UseMoles<TypeTag, TTag::ChannelNCTest> {
static constexpr bool value =
true; };
107template <
class TypeTag>
121 using GlobalPosition = Dune::FieldVector<Scalar, dimWorld>;
123 using TimeLoopPtr = std::shared_ptr<CheckPointTimeLoop<Scalar>>;
125 static constexpr auto compIdx = 1;
126 static constexpr auto transportCompIdx = Indices::conti0EqIdx + compIdx;
127 static constexpr auto transportEqIdx = Indices::conti0EqIdx + compIdx;
131 : ParentType(gridGeometry), eps_(1e-6)
133 inletVelocity_ = getParam<Scalar>(
"Problem.InletVelocity");
135 deltaP_.resize(this->gridGeometry().numCellCenterDofs());
154 {
return 273.15 + 10; }
163 return NumEqVector(0.0);
180 BoundaryTypes values;
182 if(isInlet_(globalPos))
184 values.setDirichlet(Indices::velocityXIdx);
185 values.setDirichlet(Indices::velocityYIdx);
186 values.setDirichlet(transportCompIdx);
188 values.setDirichlet(Indices::temperatureIdx);
191 else if(isOutlet_(globalPos))
193 values.setDirichlet(Indices::pressureIdx);
194 values.setOutflow(transportEqIdx);
196 values.setOutflow(Indices::energyEqIdx);
202 values.setDirichlet(Indices::velocityXIdx);
203 values.setDirichlet(Indices::velocityYIdx);
204 values.setNeumann(Indices::conti0EqIdx);
205 values.setNeumann(transportEqIdx);
207 values.setNeumann(Indices::energyEqIdx);
224 if(isInlet_(globalPos))
226 if(
time() >= 10.0 || inletVelocity_ < eps_)
228 Scalar moleFracTransportedComp = 1e-3;
230 Scalar averageMolarMassPhase = moleFracTransportedComp * FluidSystem::molarMass(compIdx)
231 + (1. - moleFracTransportedComp) * FluidSystem::molarMass(1-compIdx);
232 values[transportCompIdx] = moleFracTransportedComp * FluidSystem::molarMass(compIdx)
233 / averageMolarMassPhase;
235 values[transportCompIdx] = moleFracTransportedComp;
238 values[Indices::temperatureIdx] = 293.15;
260 PrimaryVariables values;
261 values[Indices::pressureIdx] = 1.1e+5;
262 values[transportCompIdx] = 0.0;
264 values[Indices::temperatureIdx] = 283.15;
268 values[Indices::velocityXIdx] = inletVelocity_*(globalPos[1] - this->gridGeometry().bBoxMin()[1])*(this->gridGeometry().bBoxMax()[1] - globalPos[1])
269 / (0.25*(this->gridGeometry().bBoxMax()[1] - this->gridGeometry().bBoxMin()[1])*(this->gridGeometry().bBoxMax()[1] - this->gridGeometry().bBoxMin()[1]));
271 values[Indices::velocityYIdx] = 0.0;
284 template<
class Gr
idVariables,
class SolutionVector>
287 for (
const auto& element : elements(this->gridGeometry().gridView()))
289 auto fvGeometry =
localView(this->gridGeometry());
290 fvGeometry.bindElement(element);
291 for (
auto&& scv : scvs(fvGeometry))
293 auto ccDofIdx = scv.dofIndex();
295 auto elemVolVars =
localView(gridVariables.curGridVolVars());
296 elemVolVars.bindElement(element, fvGeometry, sol);
298 deltaP_[ccDofIdx] = elemVolVars[scv].pressure() - 1.1e5;
310 timeLoop_ = timeLoop;
315 return timeLoop_->time();
320 bool isInlet_(
const GlobalPosition& globalPos)
const
322 return globalPos[0] < eps_;
325 bool isOutlet_(
const GlobalPosition& globalPos)
const
327 return globalPos[0] > this->gridGeometry().bBoxMax()[0] - eps_;
331 Scalar inletVelocity_;
332 TimeLoopPtr timeLoop_;
333 std::vector<Scalar> deltaP_;
A much simpler (and thus potentially less buggy) version of pure water.
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...
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:149
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
Property whether to use moles or kg as amount unit for balance equations.
Definition: common/properties.hh:102
The component balance index that should be replaced by the total mass/mole balance.
Definition: common/properties.hh:104
Definition: common/properties.hh:169
If disabled, the volume variables are not stored (reduces memory, but is slower)
Definition: common/properties.hh:178
specifies if data on flux vars should be saved (faster, but more memory consuming)
Definition: common/properties.hh:188
The type of the fluid system to use.
Definition: common/properties.hh:223
Switch on/off caching of face variables.
Definition: common/properties.hh:303
Navier-Stokes problem base class.
Definition: dumux/freeflow/navierstokes/problem.hh:63
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
Test problem for the one-phase (Navier-)Stokes model.
Definition: test/freeflow/navierstokesnc/channel/problem.hh:109
Scalar time() const
Definition: test/freeflow/navierstokesnc/channel/problem.hh:313
BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
Specifies which kind of boundary condition should be used for which equation on a given boundary cont...
Definition: test/freeflow/navierstokesnc/channel/problem.hh:178
NumEqVector sourceAtPos(const GlobalPosition &globalPos) const
Returns the sources within the domain.
Definition: test/freeflow/navierstokesnc/channel/problem.hh:161
auto & getDeltaP() const
Definition: test/freeflow/navierstokesnc/channel/problem.hh:303
PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
Evaluates the initial value for a control volume.
Definition: test/freeflow/navierstokesnc/channel/problem.hh:258
Scalar temperature() const
Returns the temperature within the domain in [K].
Definition: test/freeflow/navierstokesnc/channel/problem.hh:153
PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
Evaluates the boundary conditions for a Dirichlet control volume.
Definition: test/freeflow/navierstokesnc/channel/problem.hh:219
void calculateDeltaP(const GridVariables &gridVariables, const SolutionVector &sol)
Adds additional VTK output data to the VTKWriter.
Definition: test/freeflow/navierstokesnc/channel/problem.hh:285
void setTimeLoop(TimeLoopPtr timeLoop)
Definition: test/freeflow/navierstokesnc/channel/problem.hh:308
ChannelNCTestProblem(std::shared_ptr< const GridGeometry > gridGeometry)
Definition: test/freeflow/navierstokesnc/channel/problem.hh:130
bool shouldWriteRestartFile() const
Definition: test/freeflow/navierstokesnc/channel/problem.hh:143
Definition: test/freeflow/navierstokesnc/channel/problem.hh:53
std::tuple< NavierStokesNC, StaggeredFreeFlowModel > InheritsFrom
Definition: test/freeflow/navierstokesnc/channel/problem.hh:53
Dune::YaspGrid< 2 > type
Definition: test/freeflow/navierstokesnc/channel/problem.hh:73
Defines a type tag and some properties for ree-flow models using the staggered scheme....
#define ENABLECACHING
Definition: test/freeflow/navierstokesnc/channel/problem.hh:29