diff --git a/dumux/porousmediumflow/2p2c/implicit/model.hh b/dumux/porousmediumflow/2p2c/implicit/model.hh
index 2866422fc7e624444207a7fda036b7c5d1f60fb6..94096e3a90020206e8104d34e7ea1d3c28631024 100644
--- a/dumux/porousmediumflow/2p2c/implicit/model.hh
+++ b/dumux/porousmediumflow/2p2c/implicit/model.hh
@@ -96,24 +96,19 @@ class TwoPTwoCModel: public GET_PROP_TYPE(TypeTag, BaseModel)
{
// the parent class needs to access the variable switch
friend typename GET_PROP_TYPE(TypeTag, BaseModel);
-
using ParentType = typename GET_PROP_TYPE(TypeTag, BaseModel);
+
+ using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
- using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
- using PrimaryVariables = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
using VolumeVariables = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
using ElementVolumeVariables = typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables);
- using SolutionVector = typename GET_PROP_TYPE(TypeTag, SolutionVector);
- using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
-
- enum {
- numPhases = GET_PROP_VALUE(TypeTag, NumPhases),
- numComponents = GET_PROP_VALUE(TypeTag, NumComponents)
- };
- using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
- enum {
+ // phase indices
+ enum
+ {
switchIdx = Indices::switchIdx,
wPhaseIdx = Indices::wPhaseIdx,
@@ -130,19 +125,65 @@ class TwoPTwoCModel: public GET_PROP_TYPE(TypeTag, BaseModel)
formulation = GET_PROP_VALUE(TypeTag, Formulation)
};
- enum {
- dim = GridView::dimension,
- dimWorld = GridView::dimensionworld
- };
-
- using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
+ static constexpr int dim = GridView::dimension;
+ static constexpr int dimWorld = GridView::dimensionworld;
using GlobalPosition = Dune::FieldVector<Scalar, dimWorld>;
- static const bool useMoles = GET_PROP_VALUE(TypeTag, UseMoles);
- enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ static constexpr int numPhases = GET_PROP_VALUE(TypeTag, NumPhases);
+ static constexpr int numComponents = GET_PROP_VALUE(TypeTag, NumComponents);
+ static constexpr bool useMoles = GET_PROP_VALUE(TypeTag, UseMoles);
+ static constexpr bool isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox);
+
enum { dofCodim = isBox ? dim : 0 };
public:
+ /*!
+ * \brief Apply the initial conditions to the model.
+ *
+ * \param problem The object representing the problem which needs to
+ * be simulated.
+ */
+ void init(Problem& problem)
+ {
+ ParentType::init(problem);
+
+ // register standardized vtk output fields
+ auto& vtkOutputModule = problem.vtkOutputModule();
+ vtkOutputModule.addSecondaryVariable("Sn", [](const VolumeVariables& v){ return v.saturation(nPhaseIdx); });
+ vtkOutputModule.addSecondaryVariable("Sw", [](const VolumeVariables& v){ return v.saturation(wPhaseIdx); });
+ vtkOutputModule.addSecondaryVariable("pn", [](const VolumeVariables& v){ return v.pressure(nPhaseIdx); });
+ vtkOutputModule.addSecondaryVariable("pw", [](const VolumeVariables& v){ return v.pressure(wPhaseIdx); });
+ vtkOutputModule.addSecondaryVariable("pc", [](const VolumeVariables& v){ return v.capillaryPressure(); });
+ vtkOutputModule.addSecondaryVariable("rhow", [](const VolumeVariables& v){ return v.density(wPhaseIdx); });
+ vtkOutputModule.addSecondaryVariable("rhon", [](const VolumeVariables& v){ return v.density(nPhaseIdx); });
+ vtkOutputModule.addSecondaryVariable("mobW", [](const VolumeVariables& v){ return v.mobility(wPhaseIdx); });
+ vtkOutputModule.addSecondaryVariable("mobN", [](const VolumeVariables& v){ return v.mobility(nPhaseIdx); });
+
+ for (int i = 0; i < numPhases; ++i)
+ for (int j = 0; j < numComponents; ++j)
+ vtkOutputModule.addSecondaryVariable("x_" + FluidSystem::phaseName(i) + "^" + FluidSystem::componentName(j),
+ [i,j](const VolumeVariables& v){ return v.moleFraction(i,j); });
+
+ for (int i = 0; i < numPhases; ++i)
+ for (int j = 0; j < numComponents; ++j)
+ vtkOutputModule.addSecondaryVariable("X_" + FluidSystem::phaseName(i) + "^" + FluidSystem::componentName(j),
+ [i,j](const VolumeVariables& v){ return v.massFraction(i,j); });
+
+ vtkOutputModule.addSecondaryVariable("porosity", [](const VolumeVariables& v){ return v.porosity(); });
+ vtkOutputModule.addSecondaryVariable("temperature", [](const VolumeVariables& v){ return v.temperature(); });
+ }
+
+ /*!
+ * \brief Adds additional VTK output data to the VTKWriter. Function is called by the output module on every write.
+ */
+ template<class VtkOutputModule>
+ void addVtkOutputFields(VtkOutputModule& outputModule) const
+ {
+ auto& phasePresence = outputModule.createScalarField("phase presence", dofCodim);
+ for (std::size_t i = 0; i < phasePresence.size(); ++i)
+ phasePresence[i] = priVarSwitch().phasePresence(i);
+ }
+
/*!
* \brief One Newton iteration was finished.
* \param uCurrent The solution after the current Newton iteration
@@ -250,154 +291,6 @@ public:
return switchFlag_;
}
- /*!
- * \brief Append all quantities of interest which can be derived
- * from the solution of the current time step to the VTK
- * writer.
- *
- * \param sol The solution vector
- * \param writer The writer for multi-file VTK datasets
- */
- template<class MultiWriter>
- void addOutputVtkFields(const SolutionVector &sol,
- MultiWriter &writer)
- {
- using ScalarField = Dune::BlockVector<Dune::FieldVector<double, 1> >;
- // using VectorField = Dune::BlockVector<Dune::FieldVector<double, dim>>;
-
- // get the number of degrees of freedom
- auto numDofs = this->numDofs();
-
- // create the required scalar fields
- ScalarField *sN = writer.allocateManagedBuffer(numDofs);
- ScalarField *sW = writer.allocateManagedBuffer(numDofs);
- ScalarField *pn = writer.allocateManagedBuffer(numDofs);
- ScalarField *pw = writer.allocateManagedBuffer(numDofs);
- ScalarField *pc = writer.allocateManagedBuffer(numDofs);
- ScalarField *rhoW = writer.allocateManagedBuffer(numDofs);
- ScalarField *rhoN = writer.allocateManagedBuffer(numDofs);
- ScalarField *mobW = writer.allocateManagedBuffer(numDofs);
- ScalarField *mobN = writer.allocateManagedBuffer(numDofs);
- ScalarField *phasePresence = writer.allocateManagedBuffer(numDofs);
- ScalarField *massFrac[numPhases][numComponents];
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- massFrac[phaseIdx][compIdx] = writer.allocateManagedBuffer(numDofs);
- ScalarField *moleFrac[numPhases][numComponents];
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- moleFrac[phaseIdx][compIdx] = writer.allocateManagedBuffer(numDofs);
- ScalarField *temperature = writer.allocateManagedBuffer(numDofs);
- ScalarField *poro = writer.allocateManagedBuffer(numDofs);
- // VectorField *velocityN = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
- // VectorField *velocityW = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
- // ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
-
- // if (velocityOutput.enableOutput()) // check if velocity output is demanded
- // {
- // // initialize velocity fields
- // for (unsigned int i = 0; i < numDofs; ++i)
- // {
- // (*velocityN)[i] = Scalar(0);
- // (*velocityW)[i] = Scalar(0);
- // }
- // }
-
- auto numElements = this->gridView_().size(0);
- ScalarField *rank = writer.allocateManagedBuffer(numElements);
-
- for (const auto& element : elements(this->gridView_(), Dune::Partitions::interior))
- {
- auto eIdxGlobal = this->problem_().elementMapper().index(element);
- (*rank)[eIdxGlobal] = this->gridView_().comm().rank();
-
- auto fvGeometry = localView(this->globalFvGeometry());
- fvGeometry.bindElement(element);
-
- auto elemVolVars = localView(this->curGlobalVolVars());
- elemVolVars.bindElement(element, fvGeometry, this->curSol());
-
- for (auto&& scv : scvs(fvGeometry))
- {
- const auto& volVars = elemVolVars[scv];
- auto dofIdxGlobal = scv.dofIndex();
-
- (*sN)[dofIdxGlobal] = volVars.saturation(nPhaseIdx);
- (*sW)[dofIdxGlobal] = volVars.saturation(wPhaseIdx);
- (*pn)[dofIdxGlobal] = volVars.pressure(nPhaseIdx);
- (*pw)[dofIdxGlobal] = volVars.pressure(wPhaseIdx);
- (*pc)[dofIdxGlobal] = volVars.capillaryPressure();
- (*rhoW)[dofIdxGlobal] = volVars.density(wPhaseIdx);
- (*rhoN)[dofIdxGlobal] = volVars.density(nPhaseIdx);
- (*mobW)[dofIdxGlobal] = volVars.mobility(wPhaseIdx);
- (*mobN)[dofIdxGlobal] = volVars.mobility(nPhaseIdx);
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- {
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- {
- (*massFrac[phaseIdx][compIdx])[dofIdxGlobal] = volVars.massFraction(phaseIdx, compIdx);
- Valgrind::CheckDefined((*massFrac[phaseIdx][compIdx])[dofIdxGlobal][0]);
- }
- }
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- {
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- {
- (*moleFrac[phaseIdx][compIdx])[dofIdxGlobal] = volVars.moleFraction(phaseIdx, compIdx);
- Valgrind::CheckDefined((*moleFrac[phaseIdx][compIdx])[dofIdxGlobal][0]);
- }
- }
- (*poro)[dofIdxGlobal] = volVars.porosity();
- (*temperature)[dofIdxGlobal] = volVars.temperature();
- (*phasePresence)[dofIdxGlobal]= priVarSwitch().phasePresence(dofIdxGlobal);
- }
-
- // // velocity output
- // velocityOutput.calculateVelocity(*velocityW, elemVolVars, fvGeometry, element, wPhaseIdx);
- // velocityOutput.calculateVelocity(*velocityN, elemVolVars, fvGeometry, element, nPhaseIdx);
-
- } // loop over elements
-
- writer.attachDofData(*sN, "Sn", isBox);
- writer.attachDofData(*sW, "Sw", isBox);
- writer.attachDofData(*pn, "pn", isBox);
- writer.attachDofData(*pw, "pw", isBox);
- writer.attachDofData(*pc, "pc", isBox);
- writer.attachDofData(*rhoW, "rhoW", isBox);
- writer.attachDofData(*rhoN, "rhoN", isBox);
- writer.attachDofData(*mobW, "mobW", isBox);
- writer.attachDofData(*mobN, "mobN", isBox);
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- {
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- {
- std::ostringstream oss;
- oss << "X_" << FluidSystem::phaseName(phaseIdx) << "^" << FluidSystem::componentName(compIdx);
- writer.attachDofData(*massFrac[phaseIdx][compIdx], oss.str(), isBox);
- }
- }
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- {
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- {
- std::ostringstream oss;
- oss << "x_" << FluidSystem::phaseName(phaseIdx) << "^" << FluidSystem::componentName(compIdx);
- writer.attachDofData(*moleFrac[phaseIdx][compIdx], oss.str(), isBox);
- }
- }
- writer.attachDofData(*poro, "porosity", isBox);
- writer.attachDofData(*temperature, "temperature", isBox);
- writer.attachDofData(*phasePresence, "phase presence", isBox);
-
- // if (velocityOutput.enableOutput()) // check if velocity output is demanded
- // {
- // writer.attachDofData(*velocityW, "velocityW", isBox, dim);
- // writer.attachDofData(*velocityN, "velocityN", isBox, dim);
- // }
-
- writer.attachCellData(*rank, "process rank");
- }
-
/*!
* \brief Write the current solution to a restart file.
*
diff --git a/dumux/porousmediumflow/2p2c/implicit/volumevariables.hh b/dumux/porousmediumflow/2p2c/implicit/volumevariables.hh
index 801cfac4c7f3753bb511fbdb5746d637acc73be8..d116a8f68176cd3d5ae7943d4c8a3725b4d56192 100644
--- a/dumux/porousmediumflow/2p2c/implicit/volumevariables.hh
+++ b/dumux/porousmediumflow/2p2c/implicit/volumevariables.hh
@@ -46,18 +46,20 @@ class TwoPTwoCVolumeVariables : public ImplicitVolumeVariables<TypeTag>
{
using ParentType = ImplicitVolumeVariables<TypeTag>;
using Implementation = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
+
+ using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
- using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
- using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
- using PrimaryVariables = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
+ using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using MaterialLaw = typename GET_PROP_TYPE(TypeTag, MaterialLaw);
- enum {
- numPhases = GET_PROP_VALUE(TypeTag, NumPhases),
- numComponents = GET_PROP_VALUE(TypeTag, NumComponents)
- };
+ using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
+ using SpatialParams = typename GET_PROP_TYPE(TypeTag, SpatialParams);
+ using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
+ using ElementSolution = typename GET_PROP_TYPE(TypeTag, ElementSolutionVector);
- using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
- enum {
+ // component indices
+ enum
+ {
wCompIdx = Indices::wCompIdx,
nCompIdx = Indices::nCompIdx,
wPhaseIdx = Indices::wPhaseIdx,
@@ -65,41 +67,44 @@ class TwoPTwoCVolumeVariables : public ImplicitVolumeVariables<TypeTag>
};
// present phases
- enum {
+ enum
+ {
wPhaseOnly = Indices::wPhaseOnly,
nPhaseOnly = Indices::nPhaseOnly,
bothPhases = Indices::bothPhases
};
// formulations
- enum {
+ enum
+ {
formulation = GET_PROP_VALUE(TypeTag, Formulation),
pwsn = TwoPTwoCFormulation::pwsn,
pnsw = TwoPTwoCFormulation::pnsw
};
// primary variable indices
- enum {
+ enum
+ {
switchIdx = Indices::switchIdx,
pressureIdx = Indices::pressureIdx
};
- using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using Element = typename GridView::template Codim<0>::Entity;
- enum { dim = GridView::dimension};
-
- using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
- using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
+ using PermeabilityType = typename SpatialParams::PermeabilityType;
+ using ComputeFromReferencePhase = Dumux::ComputeFromReferencePhase<Scalar, FluidSystem>;
using MiscibleMultiPhaseComposition = Dumux::MiscibleMultiPhaseComposition<Scalar, FluidSystem>;
- static const bool useMoles = GET_PROP_VALUE(TypeTag, UseMoles);
- static const bool useConstraintSolver = GET_PROP_VALUE(TypeTag, UseConstraintSolver);
+
+ static constexpr bool isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox);
+ static constexpr bool useMoles = GET_PROP_VALUE(TypeTag, UseMoles);
+ static constexpr bool useKelvinEquation = GET_PROP_VALUE(TypeTag, UseKelvinEquation);
+ static constexpr bool useConstraintSolver = GET_PROP_VALUE(TypeTag, UseConstraintSolver);
static_assert(useMoles || (!useMoles && useConstraintSolver),
"if UseMoles is set false, UseConstraintSolver has to be set to true");
- static const bool useKelvinEquation = GET_PROP_VALUE(TypeTag, UseKelvinEquation);
- using ComputeFromReferencePhase = Dumux::ComputeFromReferencePhase<Scalar, FluidSystem>;
- enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
- enum { dofCodim = isBox ? dim : 0 };
+ static constexpr int dim = GridView::dimension;
+ static constexpr int numPhases = GET_PROP_VALUE(TypeTag, NumPhases);
+ static constexpr int numComponents = GET_PROP_VALUE(TypeTag, NumComponents);
+ static constexpr int dofCodim = isBox ? dim : 0;
public:
@@ -108,19 +113,19 @@ public:
/*!
* \copydoc ImplicitVolumeVariables::update
*/
- void update(const PrimaryVariables &priVars,
+ void update(const ElementSolution &elemSol,
const Problem &problem,
const Element &element,
const SubControlVolume& scv)
{
- ParentType::update(priVars, problem, element, scv);
+ ParentType::update(elemSol, problem, element, scv);
- completeFluidState(priVars, problem, element, scv, fluidState_);
+ completeFluidState(elemSol, problem, element, scv, fluidState_);
/////////////
// calculate the remaining quantities
/////////////
- const auto& materialParams = problem.spatialParams().materialLawParams(element, scv);
+ const auto& materialParams = problem.spatialParams().materialLawParams(element, scv, elemSol);
// Second instance of a parameter cache.
// Could be avoided if diffusion coefficients also
@@ -147,24 +152,26 @@ public:
nCompIdx);
}
- // porosity
- porosity_ = problem.spatialParams().porosity(scv);
+ // porosity & permeabilty
+ porosity_ = problem.spatialParams().porosity(element, scv, elemSol);
+ permeability_ = problem.spatialParams().permeability(element, scv, elemSol);
}
/*!
* \copydoc ImplicitModel::completeFluidState
* \param isOldSol Specifies whether this is the previous solution or the current one
*/
- static void completeFluidState(const PrimaryVariables& priVars,
+ static void completeFluidState(const ElementSolution& elemSol,
const Problem& problem,
const Element& element,
const SubControlVolume& scv,
FluidState& fluidState)
{
- Scalar t = ParentType::temperature(priVars, problem, element, scv);
+ Scalar t = ParentType::temperature(elemSol, problem, element, scv);
fluidState.setTemperature(t);
auto phasePresence = problem.model().priVarSwitch().phasePresence(scv.dofIndex());
+ const auto& priVars = ParentType::extractDofPriVars(elemSol, scv);
/////////////
// set the saturations
@@ -191,7 +198,7 @@ public:
/////////////
// calculate capillary pressure
- const auto& materialParams = problem.spatialParams().materialLawParams(element, scv);
+ const auto& materialParams = problem.spatialParams().materialLawParams(element, scv, elemSol);
Scalar pc = MaterialLaw::pc(materialParams, 1 - sn);
if (formulation == pwsn) {
@@ -549,6 +556,12 @@ public:
Scalar porosity() const
{ return porosity_; }
+ /*!
+ * \brief Returns the average permeability within the control volume in \f$[m^2]\f$.
+ */
+ PermeabilityType permeability() const
+ { return permeability_; }
+
/*!
* \brief Returns the binary diffusion coefficients for a phase in \f$[m^2/s]\f$.
*/
@@ -559,6 +572,7 @@ public:
protected:
Scalar porosity_; //!< Effective porosity within the control volume
+ PermeabilityType permeability_; //!< Effective permeability within the control volume
Scalar relativePermeability_[numPhases]; //!< Relative permeability within the control volume
Scalar diffCoeff_[numPhases]; //!< Binary diffusion coefficients for the phases
FluidState fluidState_;
diff --git a/test/porousmediumflow/2p2c/implicit/injectionspatialparams.hh b/test/porousmediumflow/2p2c/implicit/injectionspatialparams.hh
index 48d2209b9d7bb53d8ef5d3a14f4fa5c6c79a13ba..9d5b8918ae317cd3b2d1bdd8940c80fde81d7926 100644
--- a/test/porousmediumflow/2p2c/implicit/injectionspatialparams.hh
+++ b/test/porousmediumflow/2p2c/implicit/injectionspatialparams.hh
@@ -66,26 +66,24 @@ class InjectionSpatialParams : public ImplicitSpatialParams<TypeTag>
{
using ParentType = ImplicitSpatialParams<TypeTag>;
using Grid = typename GET_PROP_TYPE(TypeTag, Grid);
- using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
- using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
- using CoordScalar = typename Grid::ctype;
+ using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+ using MaterialLaw = typename GET_PROP_TYPE(TypeTag, MaterialLaw);
+ using MaterialLawParams = typename MaterialLaw::Params;
using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
+ using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
+ using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
- enum {
- dim=GridView::dimension,
- dimWorld=GridView::dimensionworld
- };
+ static constexpr int dim = GridView::dimension;
+ static constexpr int dimWorld = GridView::dimensionworld;
+ using CoordScalar = typename Grid::ctype;
using GlobalPosition = Dune::FieldVector<CoordScalar,dimWorld>;
- using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
- using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
- using VolumeVariables = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
using Element = typename GridView::template Codim<0>::Entity;
public:
- using MaterialLaw = typename GET_PROP_TYPE(TypeTag, MaterialLaw);
- using MaterialLawParams = typename MaterialLaw::Params;
+ using PermeabilityType = Scalar;
/*!
* \brief The constructor
@@ -124,14 +122,10 @@ public:
/*!
* \brief Returns the intrinsic permeability tensor \f$[m^2]\f$
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry of the element
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global position
*/
- Scalar intrinsicPermeability(const SubControlVolume& scv,
- const VolumeVariables& volVars) const
+ Scalar permeabilityAtPos(const GlobalPosition& globalPos) const
{
- const auto& globalPos = scv.center();
if (isFineMaterial_(globalPos))
return fineK_;
return coarseK_;
@@ -140,13 +134,10 @@ public:
/*!
* \brief Returns the porosity \f$[-]\f$
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry of the element
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global position
*/
- Scalar porosity(const SubControlVolume& scv) const
+ Scalar porosityAtPos(const GlobalPosition& globalPos) const
{
- const auto& globalPos = scv.center();
if (isFineMaterial_(globalPos))
return finePorosity_;
return coarsePorosity_;
@@ -157,14 +148,10 @@ public:
* \brief Returns the parameter object for the capillary-pressure/
* saturation material law
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry of the element
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global position
*/
- const MaterialLawParams& materialLawParams(const Element &element,
- const SubControlVolume& scv) const
+ const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
{
- const auto& globalPos = scv.center();
if (isFineMaterial_(globalPos))
return fineMaterialParams_;
return coarseMaterialParams_;
@@ -175,45 +162,30 @@ public:
*
* This is only required for non-isothermal models.
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global position
*/
- Scalar solidHeatCapacity(const Element &element,
- const SubControlVolume& scv) const
- {
- return 790; // specific heat capacity of granite [J / (kg K)]
- }
+ Scalar solidHeatCapacityAtPos(const GlobalPosition& globalPos) const
+ { return 790; /*specific heat capacity of granite [J / (kg K)]*/ }
/*!
* \brief Returns the mass density \f$[kg / m^3]\f$ of the rock matrix.
*
* This is only required for non-isothermal models.
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global position
*/
- Scalar solidDensity(const Element &element,
- const SubControlVolume& scv) const
- {
- return 2700; // density of granite [kg/m^3]
- }
+ Scalar solidDensityAtPos(const GlobalPosition& globalPos) const
+ { return 2700; /*density of granite [kg/m^3]*/ }
/*!
* \brief Returns the thermal conductivity \f$\mathrm{[W/(m K)]}\f$ of the solid
*
* This is only required for non-isothermal models.
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry of the element
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global position
*/
- Scalar solidThermalConductivity(const Element &element,
- const SubControlVolume& scv) const
- {
- return lambdaSolid_;
- }
+ Scalar solidThermalConductivityAtPos(const GlobalPosition& globalPos) const
+ { return lambdaSolid_; }
private:
bool isFineMaterial_(const GlobalPosition &globalPos) const
diff --git a/test/porousmediumflow/2p2c/implicit/waterairspatialparams.hh b/test/porousmediumflow/2p2c/implicit/waterairspatialparams.hh
index e2841d7c8a11a8587f0f8617fb1c79877f76b329..0a01a96c3af1230913126dc72fc55a5808a010e4 100644
--- a/test/porousmediumflow/2p2c/implicit/waterairspatialparams.hh
+++ b/test/porousmediumflow/2p2c/implicit/waterairspatialparams.hh
@@ -64,25 +64,19 @@ template<class TypeTag>
class WaterAirSpatialParams : public ImplicitSpatialParams<TypeTag>
{
using ParentType = ImplicitSpatialParams<TypeTag>;
- using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
- using Grid = typename GET_PROP_TYPE(TypeTag, Grid);
- using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
- using CoordScalar = typename Grid::ctype;
- enum {
- dim=GridView::dimension,
- dimWorld=GridView::dimensionworld
- };
+ using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+ using MaterialLaw = typename GET_PROP_TYPE(TypeTag, MaterialLaw);
+ using MaterialLawParams = typename MaterialLaw::Params;
+ using CoordScalar = typename GridView::ctype;
+ static constexpr int dimWorld = GridView::dimensionworld;
using GlobalPosition = Dune::FieldVector<CoordScalar,dimWorld>;
- using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
- using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
- using VolumeVariables = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
- using Element = typename GridView::template Codim<0>::Entity;
public:
- using MaterialLaw = typename GET_PROP_TYPE(TypeTag, MaterialLaw);
- using MaterialLawParams = typename MaterialLaw::Params;
+ using PermeabilityType = Scalar;
/*!
* \brief The constructor
@@ -158,14 +152,10 @@ public:
* \brief Apply the intrinsic permeability tensor to a pressure
* potential gradient.
*
- * \param element The current finite element
- * \param fvGeometry The current finite volume geometry of the element
- * \param scvIdx The index of the sub-control volume
+ * \param globalPos The global position
*/
- Scalar intrinsicPermeability(const SubControlVolume& scv,
- const VolumeVariables& volVars) const
+ Scalar permeabilityAtPos(const GlobalPosition& globalPos) const
{
- const auto& globalPos = scv.center();
if (isFineMaterial_(globalPos))
return fineK_;
return coarseK_;
@@ -174,14 +164,10 @@ public:
/*!
* \brief Define the porosity \f$[-]\f$ of the spatial parameters
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry
- * \param scvIdx The local index of the sub-control volume where
- * the porosity needs to be defined
+ * \param globalPos The global position
*/
- Scalar porosity(const SubControlVolume& scv) const
+ Scalar porosityAtPos(const GlobalPosition& globalPos) const
{
- const auto& globalPos = scv.center();
if (isFineMaterial_(globalPos))
return finePorosity_;
else
@@ -192,14 +178,10 @@ public:
/*!
* \brief return the parameter object for the Brooks-Corey material law which depends on the position
*
- * \param element The current finite element
- * \param fvGeometry The current finite volume geometry of the element
- * \param scvIdx The index of the sub-control volume
+ * \param globalPos The global position
*/
- const MaterialLawParams& materialLawParams(const Element &element,
- const SubControlVolume& scv) const
+ const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
{
- const auto& globalPos = scv.center();
if (isFineMaterial_(globalPos))
return fineMaterialParams_;
else
@@ -211,44 +193,28 @@ public:
*
* This is only required for non-isothermal models.
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global positio
*/
- Scalar solidHeatCapacity(const Element &element,
- const SubControlVolume& scv) const
- {
- return 790; // specific heat capacity of granite [J / (kg K)]
- }
+ Scalar solidHeatCapacityAtPos(const GlobalPosition& globalPos) const
+ { return 790; /*specific heat capacity of granite [J / (kg K)]*/ }
/*!
* \brief Returns the mass density \f$[kg / m^3]\f$ of the rock matrix.
*
* This is only required for non-isothermal models.
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry
- * \param scvIdx The local index of the sub-control volume
+ * \param globalPos The global position
*/
- Scalar solidDensity(const Element &element,
- const SubControlVolume& scv) const
- {
- return 2700; // density of granite [kg/m^3]
- }
+ Scalar solidDensityAtPos(const GlobalPosition& globalPos) const
+ { return 2700; /*density of granite [kg/m^3]*/ }
/*!
* \brief Returns the thermal conductivity \f$\mathrm{[W/(m K)]}\f$ of the porous material.
*
- * \param element The finite element
- * \param fvGeometry The finite volume geometry
- * \param scvIdx The local index of the sub-control volume where
- * the heat capacity needs to be defined
+ * \param globalPos The global position
*/
- Scalar solidThermalConductivity(const Element &element,
- const SubControlVolume& scv) const
- {
- return lambdaSolid_;
- }
+ Scalar solidThermalConductivityAtPos(const GlobalPosition& globalPos) const
+ { return lambdaSolid_; }
private:
bool isFineMaterial_(const GlobalPosition &globalPos) const