26 #ifndef EWOMS_STOKES_MODEL_HH
27 #define EWOMS_STOKES_MODEL_HH
40 #include <opm/material/fluidsystems/GasPhase.hpp>
41 #include <opm/material/fluidsystems/LiquidPhase.hpp>
42 #include <opm/material/components/NullComponent.hpp>
43 #include <opm/material/heatconduction/FluidConduction.hpp>
44 #include <opm/material/fluidsystems/SinglePhaseFluidSystem.hpp>
45 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
49 #include <dune/common/fvector.hh>
50 #include <dune/istl/bvector.hh>
56 template <
class TypeTag>
59 namespace Properties {
75 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
78 static const int value = GridView::dimensionworld
79 + FluidSystem::numComponents;
104 #warning "No SuperLU installed. SuperLU is the recommended linear " \
105 "solver for the Stokes models."
128 typedef Opm::FluidSystems::SinglePhase<Scalar, Fluid> type;
138 typedef Opm::LiquidPhase<Scalar, Opm::NullComponent<Scalar> > type;
147 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
150 typedef Opm::CompositionalFluidState<Scalar, FluidSystem> type;
156 Opm::FluidHeatConduction<
typename GET_PROP_TYPE(TypeTag, FluidSystem),
163 HeatConductionLawParams,
227 template <
class TypeTag>
228 class StokesModel :
public GET_PROP_TYPE(TypeTag, Discretization)
230 typedef typename GET_PROP_TYPE(TypeTag, Discretization) ParentType;
231 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
232 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
233 typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
234 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
235 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
236 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
238 enum { dimWorld = GridView::dimensionworld };
240 enum { numComponents = FluidSystem::numComponents };
242 static const int vtkFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
245 typedef typename GridView::template Codim<0>::Iterator ElementIterator;
246 typedef typename GridView::template Codim<0>::Entity Element;
250 : ParentType(simulator)
259 {
return phaseIdxQueried == phaseIdx; }
266 std::ostringstream oss;
267 if (pvIdx == Indices::pressureIdx)
269 else if (Indices::moleFrac1Idx <= pvIdx && pvIdx < Indices::moleFrac1Idx + numComponents - 1)
270 oss <<
"moleFraction^"
271 << FluidSystem::componentName(pvIdx - Indices::moleFrac1Idx + 1);
272 else if (Indices::velocity0Idx <= pvIdx && pvIdx < Indices::velocity0Idx + dimWorld)
273 oss <<
"velocity_" << pvIdx - Indices::velocity0Idx;
285 std::ostringstream oss;
286 if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx + numComponents)
288 << FluidSystem::componentName(eqIdx - Indices::conti0EqIdx);
289 else if (Indices::momentum0EqIdx <= eqIdx &&
290 eqIdx < Indices::momentum0EqIdx + dimWorld)
291 oss <<
"momentum_" << eqIdx - Indices::momentum0EqIdx;
306 if (Indices::pressureIdx == pvIdx) {
324 VtkMultiWriter *vtkWriter =
dynamic_cast<VtkMultiWriter*
>(&writer);
332 unsigned numVertices = this->gridView_.size(dimWorld);
339 ScalarBuffer *moleFraction[numComponents];
340 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
344 ElementContext elemCtx(this->simulator_);
346 ElementIterator elemIt = this->gridView().template begin<0>();
347 ElementIterator elemEndIt = this->gridView().template end<0>();
348 for (; elemIt != elemEndIt; ++elemIt) {
349 const Element& elem = *elemIt;
350 if (elem.partitionType() != Dune::InteriorEntity)
353 elemCtx.updateAll(elem);
355 int numScv = elemCtx.numPrimaryDof(0);
356 for (
int dofIdx = 0; dofIdx < numScv; ++dofIdx)
358 int globalIdx = elemCtx.globalSpaceIndex(dofIdx, 0);
359 const auto &intQuants = elemCtx.intensiveQuantities(dofIdx, 0);
360 const auto &fluidState = intQuants.fluidState();
362 pressure[globalIdx] = fluidState.pressure(phaseIdx);
363 density[globalIdx] = fluidState.density(phaseIdx);
364 temperature[globalIdx] = fluidState.temperature(phaseIdx);
365 viscosity[globalIdx] = fluidState.viscosity(phaseIdx);
372 velocity[globalIdx] = 0;
373 velocity[globalIdx] += intQuants.velocityCenter();
375 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
376 (*moleFraction[compIdx])[globalIdx] = fluidState.moleFraction(phaseIdx, compIdx);
380 std::ostringstream tmp;
383 tmp <<
"pressure_" << FluidSystem::phaseName(phaseIdx);
387 tmp <<
"density_" << FluidSystem::phaseName(phaseIdx);
390 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
392 tmp <<
"moleFraction_" << FluidSystem::phaseName(phaseIdx) <<
"^"
393 << FluidSystem::componentName(compIdx);
398 tmp <<
"temperature_" << FluidSystem::phaseName(phaseIdx);
402 tmp <<
"viscosity_" << FluidSystem::phaseName(phaseIdx);
406 tmp <<
"velocity_" << FluidSystem::phaseName(phaseIdx);
Contains the intensive quantities of the Stokes model.
Definition: stokesintensivequantities.hh:49
Declares the properties required by the Stokes model.
VectorBuffer * allocateManagedVectorBuffer(int numOuter, int numInner)
Allocate a managed buffer for a vector field.
Definition: vtkmultiwriter.hh:160
void attachVectorVertexData(VectorBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:245
bool phaseIsConsidered(int phaseIdxQueried) const
Returns true iff a fluid phase is used by the model.
Definition: stokesmodel.hh:258
std::vector< Vector > VectorBuffer
Definition: baseoutputwriter.hh:48
SET_SCALAR_PROP(NumericModel, EndTime,-1e100)
The default value for the simulation's end time.
SET_TAG_PROP(FvBaseDiscretization, LinearSolverSplice, ParallelIterativeLinearSolver)
use a parallel iterative linear solver by default
std::vector< Scalar > ScalarBuffer
Definition: baseoutputwriter.hh:47
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
void prepareOutputFields() const
Prepare the quantities relevant for the current solution to be appended to the output writers...
Definition: stokesmodel.hh:316
The local residual function for problems using the Stokes model.
Definition: stokeslocalresidual.hh:47
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:485
void attachScalarVertexData(ScalarBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:186
Base class for all problems which use the Stokes model.
SET_PROP(NumericModel, ParameterTree)
Set the ParameterTree property.
Definition: basicproperties.hh:117
Implements a boundary vector for the fully implicit (Navier-)Stokes model.
SET_INT_PROP(NumericModel, GridGlobalRefinements, 0)
The base class for all output writers.
Definition: baseoutputwriter.hh:41
SET_TYPE_PROP(NumericModel, Scalar, double)
Set the default type of scalar values to double.
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:73
Definition: baseauxiliarymodule.hh:35
ScalarBuffer * allocateManagedScalarBuffer(int numEntities)
Allocate a managed buffer for a scalar field.
Definition: vtkmultiwriter.hh:147
Contains the data which is required to calculate the mass and momentum fluxes over the face of a sub-...
Definition: stokesextensivequantities.hh:52
The base class for the vertex centered finite volume discretization scheme.
Definition: vcfvdiscretization.hh:42
Contains the intensive quantities of the Stokes model.
Contains the data which is required to calculate the mass and momentum fluxes over the face of a sub-...
std::string primaryVarName(int pvIdx) const
Given an primary variable index, return a human readable name.
Definition: stokesmodel.hh:264
The local residual function for problems using the Stokes model.
Base class for all problems which use the Stokes model.
Definition: stokesproblem.hh:46
Implements a boundary vector for the fully implicit (Navier-)Stokes model.
Definition: stokesboundaryratevector.hh:46
The primary variable and equation indices of the (Navier-)Stokes model.
Definition: stokesindices.hh:42
Scalar primaryVarWeight(int globalDofIdx, int pvIdx) const
Returns the relative weight of a primary variable for calculating relative errors.
Definition: stokesmodel.hh:301
A model for the Navier-Stokes equations.
void appendOutputFields(BaseOutputWriter &writer) const
Append the quantities relevant for the current solution to an output writer.
Definition: stokesmodel.hh:322
std::string eqName(int eqIdx) const
Given an equation index, return a human readable name.
Definition: stokesmodel.hh:283
SET_BOOL_PROP(FvBaseDiscretization, EnableVtkOutput, true)
Enable the VTK output by default.
The primary variable and equation indices of the (Navier-)Stokes model.
#define INHERITS_FROM(...)
Syntactic sugar for NEW_TYPE_TAG.
Definition: propertysystem.hh:229
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:61
StokesModel(Simulator &simulator)
Definition: stokesmodel.hh:249
A model for the Navier-Stokes equations.
Definition: stokesmodel.hh:57