OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

fix most pedantic compiler warnings in the basic infrastructure

Browse Source 
i.e., using clang 3.8 to compile the test suite with the following
flags:

```
-Weverything
-Wno-documentation
-Wno-documentation-unknown-command
-Wno-c++98-compat
-Wno-c++98-compat-pedantic
-Wno-undef
-Wno-padded
-Wno-global-constructors
-Wno-exit-time-destructors
-Wno-weak-vtables
-Wno-float-equal
```

should not produce any warnings anymore. In my opinion the only flag
which would produce beneficial warnings is -Wdocumentation. This has
not been fixed in this patch because writing documentation is left for
another day (or, more likely, year).

note that this patch consists of a heavy dose of the OPM_UNUSED macro
and plenty of static_casts (to fix signedness issues). Fixing the
singedness issues were quite a nightmare and the fact that the Dune
API is quite inconsistent in that regard was not exactly helpful. :/

Finally this patch includes quite a few formatting changes (e.g., all
occurences of 'T &t' should be changed to `T& t`) and some fixes for
minor issues which I've found during the excercise.

I've made sure that all unit tests the test suite still pass
successfully and I've made sure that flow_ebos still works for Norne
and that it did not regress w.r.t. performance.

(Note that this patch does not fix compiler warnings triggered `ebos`
and `flow_ebos` but only those caused by the basic infrastructure or
the unit tests.)

v2: fix the warnings that occur if the dune-localfunctions module is
    not available. thanks to [at]atgeirr for testing.
v3: fix dune 2.3 build issue
This commit is contained in:
Andreas Lauser
2016-11-07 15:14:07 +01:00
parent 9b53930557
commit ec4b6c82dd
20 changed files with 610 additions and 435 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
examples/problems/co2injectionflash.hh
View File

@@ -50,11 +50,11 @@ public:
* \brief Guess initial values for all quantities. * \brief Guess initial values for all quantities.
*/ */
template <class FluidState, class ComponentVector> template <class FluidState, class ComponentVector>
static void guessInitial(FluidState &fluidState, const ComponentVector &globalMolarities) static void guessInitial(FluidState& fluidState, const ComponentVector& globalMolarities)
{ {
ParentType::guessInitial(fluidState, globalMolarities); ParentType::guessInitial(fluidState, globalMolarities);
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) { for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
// pressure. use something close to the reservoir pressure as initial guess // pressure. use something close to the reservoir pressure as initial guess
fluidState.setPressure(phaseIdx, 100e5); fluidState.setPressure(phaseIdx, 100e5);
} }

72
examples/problems/co2injectionproblem.hh
View File

@@ -43,6 +43,7 @@
#include <opm/material/heatconduction/Somerton.hpp> #include <opm/material/heatconduction/Somerton.hpp>
#include <opm/material/binarycoefficients/Brine_CO2.hpp> #include <opm/material/binarycoefficients/Brine_CO2.hpp>
#include <opm/material/common/UniformTabulated2DFunction.hpp> #include <opm/material/common/UniformTabulated2DFunction.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -230,7 +231,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
Co2InjectionProblem(Simulator &simulator) Co2InjectionProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -245,11 +246,11 @@ public:
temperatureLow_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemTemperatureLow); temperatureLow_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemTemperatureLow);
temperatureHigh_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemTemperatureHigh); temperatureHigh_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemTemperatureHigh);
nTemperature_ = EWOMS_GET_PARAM(TypeTag, int, FluidSystemNumTemperature); nTemperature_ = EWOMS_GET_PARAM(TypeTag, unsigned, FluidSystemNumTemperature);
nPressure_ = EWOMS_GET_PARAM(TypeTag, int, FluidSystemNumPressure);
pressureLow_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemPressureLow); pressureLow_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemPressureLow);
pressureHigh_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemPressureHigh); pressureHigh_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemPressureHigh);
nPressure_ = EWOMS_GET_PARAM(TypeTag, unsigned, FluidSystemNumPressure);
maxDepth_ = EWOMS_GET_PARAM(TypeTag, Scalar, MaxDepth); maxDepth_ = EWOMS_GET_PARAM(TypeTag, Scalar, MaxDepth);
temperature_ = EWOMS_GET_PARAM(TypeTag, Scalar, Temperature); temperature_ = EWOMS_GET_PARAM(TypeTag, Scalar, Temperature);
@@ -374,9 +375,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (inHighTemperatureRegion_(pos)) if (inHighTemperatureRegion_(pos))
return temperature_ + 100; return temperature_ + 100;
return temperature_; return temperature_;
@@ -386,10 +387,10 @@ public:
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineK_; return fineK_;
return coarseK_; return coarseK_;
@@ -399,9 +400,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return finePorosity_; return finePorosity_;
return coarsePorosity_; return coarsePorosity_;
@@ -411,10 +412,10 @@ public:
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineMaterialParams_; return fineMaterialParams_;
return coarseMaterialParams_; return coarseMaterialParams_;
@@ -426,8 +427,9 @@ public:
* In this case, we assume the rock-matrix to be granite. * In this case, we assume the rock-matrix to be granite.
*/ */
template <class Context> template <class Context>
Scalar heatCapacitySolid(const Context &context, unsigned spaceIdx, Scalar heatCapacitySolid(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
return 790 // specific heat capacity of granite [J / (kg K)] return 790 // specific heat capacity of granite [J / (kg K)]
* 2700; // density of granite [kg/m^3] * 2700; // density of granite [kg/m^3]
@@ -438,9 +440,9 @@ public:
*/ */
template <class Context> template <class Context>
const HeatConductionLawParams & const HeatConductionLawParams &
heatConductionParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const heatConductionParams(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineHeatCondParams_; return fineHeatCondParams_;
return coarseHeatCondParams_; return coarseHeatCondParams_;
@@ -457,10 +459,10 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(pos)) { if (onLeftBoundary_(pos)) {
Opm::CompositionalFluidState<Scalar, FluidSystem> fs; Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
initialFluidState_(fs, context, spaceIdx, timeIdx); initialFluidState_(fs, context, spaceIdx, timeIdx);
@@ -505,13 +507,13 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
Opm::CompositionalFluidState<Scalar, FluidSystem> fs; Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
initialFluidState_(fs, context, spaceIdx, timeIdx); initialFluidState_(fs, context, spaceIdx, timeIdx);
// const auto &matParams = this->materialLawParams(context, spaceIdx, // const auto& matParams = this->materialLawParams(context, spaceIdx,
// timeIdx); // timeIdx);
// values.assignMassConservative(fs, matParams, /*inEquilibrium=*/true); // values.assignMassConservative(fs, matParams, /*inEquilibrium=*/true);
values.assignNaive(fs); values.assignNaive(fs);
@@ -524,18 +526,22 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
template <class Context, class FluidState> template <class Context, class FluidState>
void initialFluidState_(FluidState &fs, const Context &context, void initialFluidState_(FluidState& fs,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
////// //////
// set temperature // set temperature
@@ -556,7 +562,7 @@ private:
Scalar pl = 1e5 - densityL * this->gravity()[dim - 1] * depth; Scalar pl = 1e5 - densityL * this->gravity()[dim - 1] * depth;
Scalar pC[numPhases]; Scalar pC[numPhases];
const auto &matParams = this->materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
MaterialLaw::capillaryPressures(pC, matParams, fs); MaterialLaw::capillaryPressures(pC, matParams, fs);
fs.setPressure(liquidPhaseIdx, pl + (pC[liquidPhaseIdx] - pC[liquidPhaseIdx])); fs.setPressure(liquidPhaseIdx, pl + (pC[liquidPhaseIdx] - pC[liquidPhaseIdx]));
@@ -577,19 +583,19 @@ private:
/*setEnthalpy=*/true); /*setEnthalpy=*/true);
} }
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < eps_; } { return pos[0] < eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onInlet_(const GlobalPosition &pos) const bool onInlet_(const GlobalPosition& pos) const
{ return onRightBoundary_(pos) && (5 < pos[1]) && (pos[1] < 15); } { return onRightBoundary_(pos) && (5 < pos[1]) && (pos[1] < 15); }
bool inHighTemperatureRegion_(const GlobalPosition &pos) const bool inHighTemperatureRegion_(const GlobalPosition& pos) const
{ return (pos[0] > 20) && (pos[0] < 30) && (pos[1] > 5) && (pos[1] < 35); } { return (pos[0] > 20) && (pos[0] < 30) && (pos[1] > 5) && (pos[1] < 35); }
void computeHeatCondParams_(HeatConductionLawParams &params, Scalar poro) void computeHeatCondParams_(HeatConductionLawParams& params, Scalar poro)
{ {
Scalar lambdaWater = 0.6; Scalar lambdaWater = 0.6;
Scalar lambdaGranite = 2.8; Scalar lambdaGranite = 2.8;
@@ -603,7 +609,7 @@ private:
params.setVacuumLambda(lambdaDry); params.setVacuumLambda(lambdaDry);
} }
bool isFineMaterial_(const GlobalPosition &pos) const bool isFineMaterial_(const GlobalPosition& pos) const
{ return pos[dim - 1] > fineLayerBottom_; } { return pos[dim - 1] > fineLayerBottom_; }
DimMatrix fineK_; DimMatrix fineK_;
@@ -623,8 +629,8 @@ private:
Scalar maxDepth_; Scalar maxDepth_;
Scalar eps_; Scalar eps_;
int nTemperature_; unsigned nTemperature_;
int nPressure_; unsigned nPressure_;
Scalar pressureLow_, pressureHigh_; Scalar pressureLow_, pressureHigh_;
Scalar temperatureLow_, temperatureHigh_; Scalar temperatureLow_, temperatureHigh_;

69
examples/problems/cuvetteproblem.hh
View File

@@ -28,6 +28,8 @@
#ifndef EWOMS_CUVETTE_PROBLEM_HH #ifndef EWOMS_CUVETTE_PROBLEM_HH
#define EWOMS_CUVETTE_PROBLEM_HH #define EWOMS_CUVETTE_PROBLEM_HH
#include <ewoms/models/pvs/pvsproperties.hh>
#include <opm/material/fluidstates/CompositionalFluidState.hpp> #include <opm/material/fluidstates/CompositionalFluidState.hpp>
#include <opm/material/fluidstates/ImmiscibleFluidState.hpp> #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
#include <opm/material/fluidsystems/H2OAirMesityleneFluidSystem.hpp> #include <opm/material/fluidsystems/H2OAirMesityleneFluidSystem.hpp>
@@ -36,8 +38,8 @@
#include <opm/material/heatconduction/Somerton.hpp> #include <opm/material/heatconduction/Somerton.hpp>
#include <opm/material/constraintsolvers/MiscibleMultiPhaseComposition.hpp> #include <opm/material/constraintsolvers/MiscibleMultiPhaseComposition.hpp>
#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp> #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <ewoms/models/pvs/pvsproperties.hh> #include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -186,7 +188,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
CuvetteProblem(Simulator &simulator) CuvetteProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
, eps_(1e-6) , eps_(1e-6)
{ } { }
@@ -317,17 +319,19 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 293.15; /* [K] */ } { return 293.15; /* [K] */ }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineK_; return fineK_;
return coarseK_; return coarseK_;
@@ -337,9 +341,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return finePorosity_; return finePorosity_;
else else
@@ -350,10 +354,10 @@ public:
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineMaterialParams_; return fineMaterialParams_;
else else
@@ -365,15 +369,18 @@ public:
*/ */
template <class Context> template <class Context>
const HeatConductionLawParams & const HeatConductionLawParams &
heatConductionParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const heatConductionParams(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return heatCondParams_; } { return heatCondParams_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::heatCapacitySolid * \copydoc FvBaseMultiPhaseProblem::heatCapacitySolid
*/ */
template <class Context> template <class Context>
Scalar heatCapacitySolid(const Context &context, unsigned spaceIdx, Scalar heatCapacitySolid(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
return 850 // specific heat capacity [J / (kg K)] return 850 // specific heat capacity [J / (kg K)]
* 2650; // density of sand [kg/m^3] * 2650; // density of sand [kg/m^3]
@@ -390,10 +397,10 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onRightBoundary_(pos)) { if (onRightBoundary_(pos)) {
Opm::CompositionalFluidState<Scalar, FluidSystem> fs; Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
@@ -439,14 +446,14 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
Opm::CompositionalFluidState<Scalar, FluidSystem> fs; Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
initialFluidState_(fs, context, spaceIdx, timeIdx); initialFluidState_(fs, context, spaceIdx, timeIdx);
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
values.assignMassConservative(fs, matParams, /*inEquilibrium=*/false); values.assignMassConservative(fs, matParams, /*inEquilibrium=*/false);
} }
@@ -457,32 +464,34 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < eps_; } { return pos[0] < eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < eps_; } { return pos[1] < eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool isContaminated_(const GlobalPosition &pos) const bool isContaminated_(const GlobalPosition& pos) const
{ {
return (0.20 <= pos[0]) && (pos[0] <= 0.80) && (0.4 <= pos[1]) return (0.20 <= pos[0]) && (pos[0] <= 0.80) && (0.4 <= pos[1])
&& (pos[1] <= 0.65); && (pos[1] <= 0.65);
} }
bool isFineMaterial_(const GlobalPosition &pos) const bool isFineMaterial_(const GlobalPosition& pos) const
{ {
if (0.13 <= pos[0] && 1.20 >= pos[0] && 0.32 <= pos[1] && pos[1] <= 0.57) if (0.13 <= pos[0] && 1.20 >= pos[0] && 0.32 <= pos[1] && pos[1] <= 0.57)
return true; return true;
@@ -493,10 +502,10 @@ private:
} }
template <class FluidState, class Context> template <class FluidState, class Context>
void initialFluidState_(FluidState &fs, const Context &context, void initialFluidState_(FluidState& fs, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
fs.setTemperature(293.0 /*[K]*/); fs.setTemperature(293.0 /*[K]*/);
@@ -508,7 +517,7 @@ private:
fs.setSaturation(gasPhaseIdx, 1 - 0.12 - 0.07); fs.setSaturation(gasPhaseIdx, 1 - 0.12 - 0.07);
// set the capillary pressures // set the capillary pressures
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
Scalar pc[numPhases]; Scalar pc[numPhases];
MaterialLaw::capillaryPressures(pc, matParams, fs); MaterialLaw::capillaryPressures(pc, matParams, fs);
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
@@ -525,7 +534,7 @@ private:
fs.setSaturation(naplPhaseIdx, 0); fs.setSaturation(naplPhaseIdx, 0);
// set the capillary pressures // set the capillary pressures
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
Scalar pc[numPhases]; Scalar pc[numPhases];
MaterialLaw::capillaryPressures(pc, matParams, fs); MaterialLaw::capillaryPressures(pc, matParams, fs);
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
@@ -554,7 +563,7 @@ private:
} }
} }
void computeHeatCondParams_(HeatConductionLawParams &params, Scalar poro) void computeHeatCondParams_(HeatConductionLawParams& params, Scalar poro)
{ {
Scalar lambdaGranite = 2.8; // [W / (K m)] Scalar lambdaGranite = 2.8; // [W / (K m)]

44
examples/problems/diffusionproblem.hh
View File

@@ -30,7 +30,6 @@
#include <ewoms/models/ncp/ncpproperties.hh> #include <ewoms/models/ncp/ncpproperties.hh>
#include <dune/grid/yaspgrid.hh>
#include <ewoms/io/cubegridmanager.hh> #include <ewoms/io/cubegridmanager.hh>
#include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp> #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
@@ -38,7 +37,9 @@
#include <opm/material/fluidsystems/H2ON2FluidSystem.hpp> #include <opm/material/fluidsystems/H2ON2FluidSystem.hpp>
#include <opm/material/fluidstates/CompositionalFluidState.hpp> #include <opm/material/fluidstates/CompositionalFluidState.hpp>
#include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp> #include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh>
#include <dune/common/version.hh> #include <dune/common/version.hh>
#include <dune/common/fvector.hh> #include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh> #include <dune/common/fmatrix.hh>
@@ -174,7 +175,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
DiffusionProblem(Simulator &simulator) DiffusionProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -237,30 +238,37 @@ public:
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return K_; } { return K_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 0.35; } { return 0.35; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams&
unsigned spaceIdx, unsigned timeIdx) const materialLawParams(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return materialParams_; } { return materialParams_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } { return temperature_; }
//! \} //! \}
@@ -276,8 +284,10 @@ public:
* This problem sets no-flow boundaries everywhere. * This problem sets no-flow boundaries everywhere.
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values,
unsigned spaceIdx, unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ values.setNoFlow(); } { values.setNoFlow(); }
//! \} //! \}
@@ -291,10 +301,12 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onLeftSide_(pos)) if (onLeftSide_(pos))
values.assignNaive(leftInitialFluidState_); values.assignNaive(leftInitialFluidState_);
else else
@@ -308,14 +320,16 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLeftSide_(const GlobalPosition &pos) const bool onLeftSide_(const GlobalPosition& pos) const
{ return pos[0] < (this->boundingBoxMin()[0] + this->boundingBoxMax()[0]) / 2; } { return pos[0] < (this->boundingBoxMin()[0] + this->boundingBoxMax()[0]) / 2; }
void setupInitialFluidStates_() void setupInitialFluidStates_()

65
examples/problems/fingerproblem.hh
View File

@@ -28,9 +28,6 @@
#ifndef EWOMS_FINGER_PROBLEM_HH #ifndef EWOMS_FINGER_PROBLEM_HH
#define EWOMS_FINGER_PROBLEM_HH #define EWOMS_FINGER_PROBLEM_HH
// uncomment to run problem in 3d
// #define GRIDDIM 3
#include <ewoms/io/structuredgridmanager.hh> #include <ewoms/io/structuredgridmanager.hh>
#include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp> #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
@@ -219,7 +216,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
FingerProblem(Simulator &simulator) FingerProblem(Simulator& simulator)
: ParentType(simulator), : ParentType(simulator),
materialParams_( simulator.gridManager().grid(), codim ) materialParams_( simulator.gridManager().grid(), codim )
{ {
@@ -332,15 +329,15 @@ public:
ElementContext elemCtx(this->simulator()); ElementContext elemCtx(this->simulator());
auto elemIt = this->gridView().template begin<0>(); auto elemIt = this->gridView().template begin<0>();
const auto &elemEndIt = this->gridView().template end<0>(); const auto& elemEndIt = this->gridView().template end<0>();
for (; elemIt != elemEndIt; ++elemIt) { for (; elemIt != elemEndIt; ++elemIt) {
const auto& elem = *elemIt; const auto& elem = *elemIt;
elemCtx.updateAll( elem ); elemCtx.updateAll( elem );
const int numDofs = elemCtx.numDof(/*timeIdx=*/0); size_t numDofs = elemCtx.numDof(/*timeIdx=*/0);
for (int scvIdx = 0; scvIdx < numDofs; ++scvIdx) for (unsigned scvIdx = 0; scvIdx < numDofs; ++scvIdx)
{ {
MaterialLawParams& materialParam = materialLawParams( elemCtx, scvIdx, /*timeIdx=*/0 ); MaterialLawParams& materialParam = materialLawParams( elemCtx, scvIdx, /*timeIdx=*/0 );
const auto &fs = elemCtx.intensiveQuantities(scvIdx, /*timeIdx=*/0).fluidState(); const auto& fs = elemCtx.intensiveQuantities(scvIdx, /*timeIdx=*/0).fluidState();
ParkerLenhard::update(materialParam, fs); ParkerLenhard::update(materialParam, fs);
} }
} }
@@ -357,46 +354,45 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& /*context*/, unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
{ return temperature_; } { return temperature_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& /*context*/, unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
unsigned timeIdx) const
{ return K_; } { return K_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& /*context*/, unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
{ return 0.4; } { return 0.4; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
MaterialLawParams &materialLawParams(const Context &context, MaterialLawParams& materialLawParams(const Context& context,
const int spaceIdx, const int timeIdx) unsigned spaceIdx, unsigned timeIdx)
{ {
const auto& entity = context.stencil(timeIdx).entity( spaceIdx ); const auto& entity = context.stencil(timeIdx).entity(spaceIdx);
assert( materialParams_[ entity ] ); assert(materialParams_[entity]);
return *(materialParams_[ entity ] ); return *materialParams_[entity];
} }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& context,
const int spaceIdx, const int timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const auto& entity = context.stencil(timeIdx).entity( spaceIdx ); const auto& entity = context.stencil(timeIdx).entity( spaceIdx );
assert( materialParams_[ entity ] ); assert(materialParams_[entity]);
return *(materialParams_[ entity ] ); return *materialParams_[entity];
} }
//! \} //! \}
@@ -410,10 +406,10 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(pos) || onRightBoundary_(pos) || onLowerBoundary_(pos)) if (onLeftBoundary_(pos) || onRightBoundary_(pos) || onLowerBoundary_(pos))
values.setNoFlow(); values.setNoFlow();
@@ -441,8 +437,7 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values, const Context& /*context*/, unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
unsigned timeIdx) const
{ {
// assign the primary variables // assign the primary variables
values.assignNaive(initialFluidState_); values.assignNaive(initialFluidState_);
@@ -452,10 +447,10 @@ public:
* \copydoc FvBaseProblem::constraints * \copydoc FvBaseProblem::constraints
*/ */
template <class Context> template <class Context>
void constraints(Constraints &constraints, const Context &context, void constraints(Constraints& constraints, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (onUpperBoundary_(pos) && !onInlet_(pos)) { if (onUpperBoundary_(pos) && !onInlet_(pos)) {
constraints.setActive(true); constraints.setActive(true);
@@ -474,25 +469,25 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate, const Context& /*context*/,
unsigned timeIdx) const unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < this->boundingBoxMin()[0] + eps_; } { return pos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < this->boundingBoxMin()[1] + eps_; } { return pos[1] < this->boundingBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool onInlet_(const GlobalPosition &pos) const bool onInlet_(const GlobalPosition& pos) const
{ {
Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0]; Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0];
Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width; Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width;
@@ -512,7 +507,7 @@ private:
void setupInitialFluidState_() void setupInitialFluidState_()
{ {
auto &fs = initialFluidState_; auto& fs = initialFluidState_;
fs.setPressure(wettingPhaseIdx, /*pressure=*/1e5); fs.setPressure(wettingPhaseIdx, /*pressure=*/1e5);
Scalar Sw = EWOMS_GET_PARAM(TypeTag, Scalar, InitialWaterSaturation); Scalar Sw = EWOMS_GET_PARAM(TypeTag, Scalar, InitialWaterSaturation);

88
examples/problems/fractureproblem.hh
View File

@@ -28,7 +28,6 @@
#ifndef EWOMS_FRACTURE_PROBLEM_HH #ifndef EWOMS_FRACTURE_PROBLEM_HH
#define EWOMS_FRACTURE_PROBLEM_HH #define EWOMS_FRACTURE_PROBLEM_HH
#if HAVE_DUNE_ALUGRID #if HAVE_DUNE_ALUGRID
// avoid reordering of macro elements, otherwise this problem won't work // avoid reordering of macro elements, otherwise this problem won't work
#define DISABLE_ALUGRID_SFC_ORDERING 1 #define DISABLE_ALUGRID_SFC_ORDERING 1
@@ -40,16 +39,17 @@
#include <ewoms/models/discretefracture/discretefracturemodel.hh> #include <ewoms/models/discretefracture/discretefracturemodel.hh>
#include <ewoms/io/dgfgridmanager.hh> #include <ewoms/io/dgfgridmanager.hh>
#include <opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp> #include <opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp>
#include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp> #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
#include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp> #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
#include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp> #include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp>
#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp> #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
#include <opm/material/heatconduction/Somerton.hpp> #include <opm/material/heatconduction/Somerton.hpp>
#include <opm/material/fluidsystems/TwoPhaseImmiscibleFluidSystem.hpp> #include <opm/material/fluidsystems/TwoPhaseImmiscibleFluidSystem.hpp>
#include <opm/material/components/SimpleH2O.hpp> #include <opm/material/components/SimpleH2O.hpp>
#include <opm/material/components/Dnapl.hpp> #include <opm/material/components/Dnapl.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/common/version.hh> #include <dune/common/version.hh>
#include <dune/common/fmatrix.hh> #include <dune/common/fmatrix.hh>
@@ -220,7 +220,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
FractureProblem(Simulator &simulator) FractureProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -316,7 +316,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } { return temperature_; }
// \} // \}
@@ -330,8 +332,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return matrixK_; } { return matrixK_; }
/*! /*!
@@ -340,16 +343,18 @@ public:
* \copydoc Doxygen::contextParams * \copydoc Doxygen::contextParams
*/ */
template <class Context> template <class Context>
const DimMatrix &fractureIntrinsicPermeability(const Context &context, const DimMatrix& fractureIntrinsicPermeability(const Context& OPM_UNUSED context,
unsigned spaceIdx, unsigned OPM_UNUSED spaceIdx,
unsigned timeIdx) const unsigned OPM_UNUSED timeIdx) const
{ return fractureK_; } { return fractureK_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return matrixPorosity_; } { return matrixPorosity_; }
/*! /*!
@@ -358,16 +363,18 @@ public:
* \copydoc Doxygen::contextParams * \copydoc Doxygen::contextParams
*/ */
template <class Context> template <class Context>
Scalar fracturePorosity(const Context &context, unsigned spaceIdx, Scalar fracturePorosity(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return fracturePorosity_; } { return fracturePorosity_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& OPM_UNUSED context,
unsigned spaceIdx, unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return matrixMaterialParams_; } { return matrixMaterialParams_; }
/*! /*!
@@ -376,15 +383,15 @@ public:
* \copydoc Doxygen::contextParams * \copydoc Doxygen::contextParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &fractureMaterialLawParams(const Context &context, const MaterialLawParams& fractureMaterialLawParams(const Context& OPM_UNUSED context,
unsigned spaceIdx, unsigned OPM_UNUSED spaceIdx,
unsigned timeIdx) const unsigned OPM_UNUSED timeIdx) const
{ return fractureMaterialParams_; } { return fractureMaterialParams_; }
/*! /*!
* \brief Returns the object representating the fracture topology. * \brief Returns the object representating the fracture topology.
*/ */
const FractureMapper &fractureMapper() const const FractureMapper& fractureMapper() const
{ return this->simulator().gridManager().fractureMapper(); } { return this->simulator().gridManager().fractureMapper(); }
/*! /*!
@@ -400,8 +407,10 @@ public:
* \param timeIdx The index used by the time discretization. * \param timeIdx The index used by the time discretization.
*/ */
template <class Context> template <class Context>
Scalar fractureWidth(const Context &context, unsigned spaceIdx1, unsigned spaceIdx2, Scalar fractureWidth(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx1,
unsigned OPM_UNUSED spaceIdx2,
unsigned OPM_UNUSED timeIdx) const
{ return fractureWidth_; } { return fractureWidth_; }
/*! /*!
@@ -409,7 +418,9 @@ public:
*/ */
template <class Context> template <class Context>
const HeatConductionLawParams & const HeatConductionLawParams &
heatConductionParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const heatConductionParams(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return heatCondParams_; } { return heatCondParams_; }
/*! /*!
@@ -418,8 +429,9 @@ public:
* In this case, we assume the rock-matrix to be granite. * In this case, we assume the rock-matrix to be granite.
*/ */
template <class Context> template <class Context>
Scalar heatCapacitySolid(const Context &context, unsigned spaceIdx, Scalar heatCapacitySolid(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
return 790 // specific heat capacity of granite [J / (kg K)] return 790 // specific heat capacity of granite [J / (kg K)]
* 2700; // density of granite [kg/m^3] * 2700; // density of granite [kg/m^3]
@@ -436,10 +448,10 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (onRightBoundary_(pos)) { if (onRightBoundary_(pos)) {
// on the right boundary, we impose a free-flow // on the right boundary, we impose a free-flow
@@ -474,10 +486,10 @@ public:
* \copydoc FvBaseProblem::constraints * \copydoc FvBaseProblem::constraints
*/ */
template <class Context> template <class Context>
void constraints(Constraints &constraints, const Context &context, void constraints(Constraints& constraints, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (!onLeftBoundary_(pos)) if (!onLeftBoundary_(pos))
// only impose constraints adjacent to the left boundary // only impose constraints adjacent to the left boundary
@@ -520,8 +532,10 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
FluidState fluidState; FluidState fluidState;
fluidState.setTemperature(temperature_); fluidState.setTemperature(temperature_);
@@ -542,26 +556,28 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
// \} // \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < this->boundingBoxMin()[0] + eps_; } { return pos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < this->boundingBoxMin()[1] + eps_; } { return pos[1] < this->boundingBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
void computeHeatCondParams_(HeatConductionLawParams &params, Scalar poro) void computeHeatCondParams_(HeatConductionLawParams& params, Scalar poro)
{ {
Scalar lambdaGranite = 2.8; // [W / (K m)] Scalar lambdaGranite = 2.8; // [W / (K m)]

44
examples/problems/groundwaterproblem.hh
View File

@@ -34,6 +34,7 @@
#include <opm/material/components/SimpleH2O.hpp> #include <opm/material/components/SimpleH2O.hpp>
#include <opm/material/fluidstates/ImmiscibleFluidState.hpp> #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
#include <opm/material/fluidsystems/LiquidPhase.hpp> #include <opm/material/fluidsystems/LiquidPhase.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -157,7 +158,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
GroundWaterProblem(Simulator &simulator) GroundWaterProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -263,25 +264,32 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 273.15 + 10; } // 10C { return 273.15 + 10; } // 10C
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 0.4; } { return 0.4; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
return isInLens_(context.pos(spaceIdx, timeIdx)) ? intrinsicPermLens_ if (isInLens_(context.pos(spaceIdx, timeIdx)))
: intrinsicPerm_; return intrinsicPermLens_;
else
return intrinsicPerm_;
} }
//! \} //! \}
@@ -294,10 +302,10 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &globalPos = context.pos(spaceIdx, timeIdx); const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
if (onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos)) { if (onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos)) {
Scalar pressure; Scalar pressure;
@@ -333,10 +341,12 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
// const GlobalPosition &globalPos = context.pos(spaceIdx, timeIdx); // const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
values[pressure0Idx] = 1.0e+5; // + 9.81*1.23*(20-globalPos[dim-1]); values[pressure0Idx] = 1.0e+5; // + 9.81*1.23*(20-globalPos[dim-1]);
} }
@@ -344,20 +354,22 @@ public:
* \copydoc FvBaseProblem::source * \copydoc FvBaseProblem::source
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[dim - 1] < eps_; } { return pos[dim - 1] < eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[dim - 1] > this->boundingBoxMax()[dim - 1] - eps_; } { return pos[dim - 1] > this->boundingBoxMax()[dim - 1] - eps_; }
bool isInLens_(const GlobalPosition &pos) const bool isInLens_(const GlobalPosition& pos) const
{ {
return lensLowerLeft_[0] <= pos[0] && pos[0] <= lensUpperRight_[0] return lensLowerLeft_[0] <= pos[0] && pos[0] <= lensUpperRight_[0]
&& lensLowerLeft_[1] <= pos[1] && pos[1] <= lensUpperRight_[1]; && lensLowerLeft_[1] <= pos[1] && pos[1] <= lensUpperRight_[1];

67
examples/problems/infiltrationproblem.hh
View File

@@ -35,6 +35,8 @@
#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp> #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
#include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp> #include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp>
#include <opm/material/heatconduction/Somerton.hpp> #include <opm/material/heatconduction/Somerton.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -189,7 +191,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
InfiltrationProblem(Simulator &simulator) InfiltrationProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
, eps_(1e-6) , eps_(1e-6)
{ } { }
@@ -277,17 +279,21 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } { return temperature_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix&
unsigned timeIdx) const intrinsicPermeability(const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineK_; return fineK_;
return coarseK_; return coarseK_;
@@ -297,15 +303,19 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return porosity_; } { return porosity_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams&
unsigned spaceIdx, unsigned timeIdx) const materialLawParams(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return materialParams_; } { return materialParams_; }
/*! /*!
@@ -314,8 +324,9 @@ public:
* In this case, we assume the rock-matrix to be quartz. * In this case, we assume the rock-matrix to be quartz.
*/ */
template <class Context> template <class Context>
Scalar heatCapacitySolid(const Context &context, unsigned spaceIdx, Scalar heatCapacitySolid(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
return 850. // specific heat capacity [J / (kg K)] return 850. // specific heat capacity [J / (kg K)]
* 2650.; // density of sand [kg/m^3] * 2650.; // density of sand [kg/m^3]
@@ -332,10 +343,12 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(pos) || onRightBoundary_(pos)) { if (onLeftBoundary_(pos) || onRightBoundary_(pos)) {
Opm::CompositionalFluidState<Scalar, FluidSystem> fs; Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
@@ -366,14 +379,16 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
Opm::CompositionalFluidState<Scalar, FluidSystem> fs; Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
initialFluidState_(fs, context, spaceIdx, timeIdx); initialFluidState_(fs, context, spaceIdx, timeIdx);
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
values.assignMassConservative(fs, matParams, /*inEquilibrium=*/true); values.assignMassConservative(fs, matParams, /*inEquilibrium=*/true);
Valgrind::CheckDefined(values); Valgrind::CheckDefined(values);
} }
@@ -385,30 +400,32 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < eps_; } { return pos[0] < eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < eps_; } { return pos[1] < eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool onInlet_(const GlobalPosition &pos) const bool onInlet_(const GlobalPosition& pos) const
{ return onUpperBoundary_(pos) && 50 < pos[0] && pos[0] < 75; } { return onUpperBoundary_(pos) && 50 < pos[0] && pos[0] < 75; }
template <class FluidState, class Context> template <class FluidState, class Context>
void initialFluidState_(FluidState &fs, const Context &context, void initialFluidState_(FluidState& fs, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition pos = context.pos(spaceIdx, timeIdx); const GlobalPosition pos = context.pos(spaceIdx, timeIdx);
@@ -421,7 +438,7 @@ private:
pc = 0.0; pc = 0.0;
// set pressures // set pressures
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
Scalar Sw = matParams.Swr(); Scalar Sw = matParams.Swr();
Scalar Swr = matParams.Swr(); Scalar Swr = matParams.Swr();
Scalar Sgr = matParams.Sgr(); Scalar Sgr = matParams.Sgr();
@@ -466,7 +483,7 @@ private:
1 - fs.moleFraction(waterPhaseIdx, H2OIdx)); 1 - fs.moleFraction(waterPhaseIdx, H2OIdx));
} }
bool isFineMaterial_(const GlobalPosition &pos) const bool isFineMaterial_(const GlobalPosition& pos) const
{ return 70. <= pos[0] && pos[0] <= 85. && 7.0 <= pos[1] && pos[1] <= 7.50; } { return 70. <= pos[0] && pos[0] <= 85. && 7.0 <= pos[1] && pos[1] <= 7.50; }
DimMatrix fineK_; DimMatrix fineK_;

53
examples/problems/lensproblem.hh
View File

@@ -40,6 +40,7 @@
#include <opm/material/fluidstates/ImmiscibleFluidState.hpp> #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
#include <opm/material/components/SimpleH2O.hpp> #include <opm/material/components/SimpleH2O.hpp>
#include <opm/material/components/Dnapl.hpp> #include <opm/material/components/Dnapl.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/common/version.hh> #include <dune/common/version.hh>
#include <dune/common/fvector.hh> #include <dune/common/fvector.hh>
@@ -223,7 +224,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
LensProblem(Simulator &simulator) LensProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -311,10 +312,10 @@ public:
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &globalPos = context.pos(spaceIdx, timeIdx); const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
if (isInLens_(globalPos)) if (isInLens_(globalPos))
return lensK_; return lensK_;
@@ -325,17 +326,19 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 0.4; } { return 0.4; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &globalPos = context.pos(spaceIdx, timeIdx); const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
if (isInLens_(globalPos)) if (isInLens_(globalPos))
return lensMaterialParams_; return lensMaterialParams_;
@@ -346,7 +349,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } { return temperature_; }
//! \} //! \}
@@ -409,10 +414,12 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, void boundary(BoundaryRateVector& values,
const Context &context, unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(pos) || onRightBoundary_(pos)) { if (onLeftBoundary_(pos) || onRightBoundary_(pos)) {
// free flow boundary // free flow boundary
@@ -441,7 +448,7 @@ public:
} }
// specify a full fluid state using pw and Sw // specify a full fluid state using pw and Sw
const MaterialLawParams &matParams = this->materialLawParams(context, spaceIdx, timeIdx); const MaterialLawParams& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
Opm::ImmiscibleFluidState<Scalar, FluidSystem, Opm::ImmiscibleFluidState<Scalar, FluidSystem,
/*storeEnthalpy=*/false> fs; /*storeEnthalpy=*/false> fs;
@@ -482,9 +489,9 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
Scalar depth = this->boundingBoxMax()[1] - pos[1]; Scalar depth = this->boundingBoxMax()[1] - pos[1];
Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs; Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs;
@@ -504,7 +511,7 @@ public:
Scalar pw = 1e5 - densityW * this->gravity()[1] * depth; Scalar pw = 1e5 - densityW * this->gravity()[1] * depth;
// calculate the capillary pressure // calculate the capillary pressure
const MaterialLawParams &matParams = this->materialLawParams(context, spaceIdx, timeIdx); const MaterialLawParams& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
Scalar pC[numPhases]; Scalar pC[numPhases];
MaterialLaw::capillaryPressures(pC, matParams, fs); MaterialLaw::capillaryPressures(pC, matParams, fs);
@@ -523,14 +530,16 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool isInLens_(const GlobalPosition &pos) const bool isInLens_(const GlobalPosition& pos) const
{ {
for (unsigned i = 0; i < dim; ++i) { for (unsigned i = 0; i < dim; ++i) {
if (pos[i] < lensLowerLeft_[i] - eps_ || pos[i] > lensUpperRight_[i] if (pos[i] < lensLowerLeft_[i] - eps_ || pos[i] > lensUpperRight_[i]
@@ -540,19 +549,19 @@ private:
return true; return true;
} }
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < this->boundingBoxMin()[0] + eps_; } { return pos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < this->boundingBoxMin()[1] + eps_; } { return pos[1] < this->boundingBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool onInlet_(const GlobalPosition &pos) const bool onInlet_(const GlobalPosition& pos) const
{ {
Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0]; Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0];
Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width; Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width;

72
examples/problems/obstacleproblem.hh
View File

@@ -38,6 +38,7 @@
#include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp> #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp> #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
#include <opm/material/heatconduction/Somerton.hpp> #include <opm/material/heatconduction/Somerton.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -178,7 +179,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
ObstacleProblem(Simulator &simulator) ObstacleProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -195,11 +196,11 @@ public:
// initialize the tables of the fluid system // initialize the tables of the fluid system
Scalar Tmin = temperature_ - 1.0; Scalar Tmin = temperature_ - 1.0;
Scalar Tmax = temperature_ + 1.0; Scalar Tmax = temperature_ + 1.0;
int nT = 3; unsigned nT = 3;
Scalar pmin = 1.0e5 * 0.75; Scalar pmin = 1.0e5 * 0.75;
Scalar pmax = 2.0e5 * 1.25; Scalar pmax = 2.0e5 * 1.25;
int np = 1000; unsigned np = 1000;
FluidSystem::init(Tmin, Tmax, nT, pmin, pmax, np); FluidSystem::init(Tmin, Tmax, nT, pmin, pmax, np);
@@ -298,15 +299,19 @@ public:
* This problem simply assumes a constant temperature. * This problem simply assumes a constant temperature.
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } { return temperature_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix&
unsigned timeIdx) const intrinsicPermeability(const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
if (isFineMaterial_(context.pos(spaceIdx, timeIdx))) if (isFineMaterial_(context.pos(spaceIdx, timeIdx)))
return fineK_; return fineK_;
@@ -317,9 +322,11 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return finePorosity_; return finePorosity_;
else else
@@ -330,10 +337,12 @@ public:
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams&
unsigned spaceIdx, unsigned timeIdx) const materialLawParams(const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineMaterialParams_; return fineMaterialParams_;
else else
@@ -347,8 +356,9 @@ public:
* medium is granite. * medium is granite.
*/ */
template <class Context> template <class Context>
Scalar heatCapacitySolid(const Context &context, unsigned spaceIdx, Scalar heatCapacitySolid(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
return 790 // specific heat capacity of granite [J / (kg K)] return 790 // specific heat capacity of granite [J / (kg K)]
* 2700; // density of granite [kg/m^3] * 2700; // density of granite [kg/m^3]
@@ -359,9 +369,11 @@ public:
*/ */
template <class Context> template <class Context>
const HeatConductionLawParams & const HeatConductionLawParams &
heatConductionParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const heatConductionParams(const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineHeatCondParams_; return fineHeatCondParams_;
return coarseHeatCondParams_; return coarseHeatCondParams_;
@@ -378,10 +390,12 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onInlet_(pos)) if (onInlet_(pos))
values.setFreeFlow(context, spaceIdx, timeIdx, inletFluidState_); values.setFreeFlow(context, spaceIdx, timeIdx, inletFluidState_);
@@ -402,10 +416,12 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
values.assignMassConservative(outletFluidState_, matParams); values.assignMassConservative(outletFluidState_, matParams);
} }
@@ -416,8 +432,10 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = 0.0; } { rate = 0.0; }
//! \} //! \}
@@ -427,17 +445,17 @@ private:
* \brief Returns whether a given global position is in the * \brief Returns whether a given global position is in the
* fine-permeability region or not. * fine-permeability region or not.
*/ */
bool isFineMaterial_(const GlobalPosition &pos) const bool isFineMaterial_(const GlobalPosition& pos) const
{ return 10 <= pos[0] && pos[0] <= 20 && 0 <= pos[1] && pos[1] <= 35; } { return 10 <= pos[0] && pos[0] <= 20 && 0 <= pos[1] && pos[1] <= 35; }
bool onInlet_(const GlobalPosition &globalPos) const bool onInlet_(const GlobalPosition& globalPos) const
{ {
Scalar x = globalPos[0]; Scalar x = globalPos[0];
Scalar y = globalPos[1]; Scalar y = globalPos[1];
return x >= 60 - eps_ && y <= 10; return x >= 60 - eps_ && y <= 10;
} }
bool onOutlet_(const GlobalPosition &globalPos) const bool onOutlet_(const GlobalPosition& globalPos) const
{ {
Scalar x = globalPos[0]; Scalar x = globalPos[0];
Scalar y = globalPos[1]; Scalar y = globalPos[1];
@@ -453,7 +471,7 @@ private:
} }
template <class FluidState> template <class FluidState>
void initFluidState_(FluidState &fs, const MaterialLawParams &matParams, bool isInlet) void initFluidState_(FluidState& fs, const MaterialLawParams& matParams, bool isInlet)
{ {
unsigned refPhaseIdx; unsigned refPhaseIdx;
unsigned otherPhaseIdx; unsigned otherPhaseIdx;
@@ -512,7 +530,7 @@ private:
/*setEnthalpy=*/false); /*setEnthalpy=*/false);
} }
void computeHeatCondParams_(HeatConductionLawParams &params, Scalar poro) void computeHeatCondParams_(HeatConductionLawParams& params, Scalar poro)
{ {
Scalar lambdaWater = 0.6; Scalar lambdaWater = 0.6;
Scalar lambdaGranite = 2.8; Scalar lambdaGranite = 2.8;

40
examples/problems/outflowproblem.hh
View File

@@ -31,6 +31,7 @@
#include <opm/material/fluidstates/CompositionalFluidState.hpp> #include <opm/material/fluidstates/CompositionalFluidState.hpp>
#include <opm/material/fluidsystems/H2ON2LiquidPhaseFluidSystem.hpp> #include <opm/material/fluidsystems/H2ON2LiquidPhaseFluidSystem.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -135,7 +136,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
OutflowProblem(Simulator &simulator) OutflowProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
, eps_(1e-6) , eps_(1e-6)
{ } { }
@@ -194,7 +195,9 @@ public:
* This problem assumes a temperature. * This problem assumes a temperature.
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } // in [K] { return temperature_; } // in [K]
/*! /*!
@@ -203,8 +206,9 @@ public:
* This problem uses a constant intrinsic permeability. * This problem uses a constant intrinsic permeability.
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return perm_; } { return perm_; }
/*! /*!
@@ -213,7 +217,9 @@ public:
* This problem uses a constant porosity. * This problem uses a constant porosity.
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return porosity_; } { return porosity_; }
#if 0 #if 0
@@ -222,7 +228,7 @@ public:
* *
*/ */
template <class Context> template <class Context>
Scalar tortuosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar tortuosity(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ return tortuosity_; } { return tortuosity_; }
/*! /*!
@@ -230,7 +236,7 @@ public:
* *
*/ */
template <class Context> template <class Context>
Scalar dispersivity(const Context &context, Scalar dispersivity(const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ return 0; } { return 0; }
#endif #endif
@@ -246,10 +252,10 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &globalPos = context.pos(spaceIdx, timeIdx); const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(globalPos)) { if (onLeftBoundary_(globalPos)) {
Opm::CompositionalFluidState<Scalar, FluidSystem, Opm::CompositionalFluidState<Scalar, FluidSystem,
@@ -288,7 +294,9 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
Opm::CompositionalFluidState<Scalar, FluidSystem, /*storeEnthalpy=*/false> fs; Opm::CompositionalFluidState<Scalar, FluidSystem, /*storeEnthalpy=*/false> fs;
@@ -304,21 +312,23 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < eps_; } { return pos[0] < eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
template <class FluidState, class Context> template <class FluidState, class Context>
void initialFluidState_(FluidState &fs, const Context &context, void initialFluidState_(FluidState& fs, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
Scalar T = temperature(context, spaceIdx, timeIdx); Scalar T = temperature(context, spaceIdx, timeIdx);

56
examples/problems/powerinjectionproblem.hh
View File

@@ -28,6 +28,9 @@
#ifndef EWOMS_POWER_INJECTION_PROBLEM_HH #ifndef EWOMS_POWER_INJECTION_PROBLEM_HH
#define EWOMS_POWER_INJECTION_PROBLEM_HH #define EWOMS_POWER_INJECTION_PROBLEM_HH
#include <ewoms/models/immiscible/immisciblemodel.hh>
#include <ewoms/io/cubegridmanager.hh>
#include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp> #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
#include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp> #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
#include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp> #include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp>
@@ -36,10 +39,9 @@
#include <opm/material/fluidstates/ImmiscibleFluidState.hpp> #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
#include <opm/material/components/SimpleH2O.hpp> #include <opm/material/components/SimpleH2O.hpp>
#include <opm/material/components/Air.hpp> #include <opm/material/components/Air.hpp>
#include <ewoms/models/immiscible/immisciblemodel.hh> #include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <ewoms/io/cubegridmanager.hh>
#include <dune/common/version.hh> #include <dune/common/version.hh>
#include <dune/common/fvector.hh> #include <dune/common/fvector.hh>
@@ -193,7 +195,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
PowerInjectionProblem(Simulator &simulator) PowerInjectionProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -269,38 +271,46 @@ public:
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return K_; } { return K_; }
/*! /*!
* \copydoc ForchheimerBaseProblem::ergunCoefficient * \copydoc ForchheimerBaseProblem::ergunCoefficient
*/ */
template <class Context> template <class Context>
Scalar ergunCoefficient(const Context &context, unsigned spaceIdx, Scalar ergunCoefficient(const Context& OPM_UNUSED context,
unsigned timeIdx) const unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 0.3866; } { return 0.3866; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 0.558; } { return 0.558; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams&
unsigned spaceIdx, unsigned timeIdx) const materialLawParams(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return materialParams_; } { return materialParams_; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } { return temperature_; }
//! \} //! \}
@@ -317,10 +327,12 @@ public:
* left and a free-flow boundary on the right. * left and a free-flow boundary on the right.
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(pos)) { if (onLeftBoundary_(pos)) {
RateVector massRate(0.0); RateVector massRate(0.0);
@@ -348,8 +360,10 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
// assign the primary variables // assign the primary variables
values.assignNaive(initialFluidState_); values.assignNaive(initialFluidState_);
@@ -362,17 +376,19 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < this->boundingBoxMin()[0] + eps_; } { return pos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
void setupInitialFluidState_() void setupInitialFluidState_()

63
examples/problems/reservoirproblem.hh
View File

@@ -35,10 +35,10 @@
#include <opm/material/fluidstates/CompositionalFluidState.hpp> #include <opm/material/fluidstates/CompositionalFluidState.hpp>
#include <opm/material/fluidsystems/BlackOilFluidSystem.hpp> #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
#include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp> #include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp>
#include <opm/material/fluidsystems/blackoilpvt/DryGasPvt.hpp> #include <opm/material/fluidsystems/blackoilpvt/DryGasPvt.hpp>
#include <opm/material/fluidsystems/blackoilpvt/LiveOilPvt.hpp> #include <opm/material/fluidsystems/blackoilpvt/LiveOilPvt.hpp>
#include <opm/material/fluidsystems/blackoilpvt/ConstantCompressibilityWaterPvt.hpp> #include <opm/material/fluidsystems/blackoilpvt/ConstantCompressibilityWaterPvt.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -201,7 +201,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
ReservoirProblem(Simulator &simulator) ReservoirProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -351,10 +351,10 @@ public:
const auto& eEndIt = this->simulator().gridView().template end<0>(); const auto& eEndIt = this->simulator().gridView().template end<0>();
for (; eIt != eEndIt; ++eIt) { for (; eIt != eEndIt; ++eIt) {
elemCtx.updateStencil(*eIt); elemCtx.updateStencil(*eIt);
int nDof = elemCtx.numPrimaryDof(/*timeIdx=*/0); size_t nDof = elemCtx.numPrimaryDof(/*timeIdx=*/0);
for (int dofIdx = 0; dofIdx < nDof; ++ dofIdx) { for (unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
int globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0); unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
const GlobalPosition &pos = elemCtx.pos(dofIdx, /*timeIdx=*/0); const GlobalPosition& pos = elemCtx.pos(dofIdx, /*timeIdx=*/0);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
materialParams_[globalDofIdx] = &fineMaterialParams_; materialParams_[globalDofIdx] = &fineMaterialParams_;
@@ -431,10 +431,10 @@ public:
* above one with low permeability. * above one with low permeability.
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineK_; return fineK_;
return coarseK_; return coarseK_;
@@ -444,9 +444,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return finePorosity_; return finePorosity_;
return coarsePorosity_; return coarsePorosity_;
@@ -456,14 +456,14 @@ public:
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
unsigned globalIdx = context.globalSpaceIndex(spaceIdx, timeIdx); unsigned globalIdx = context.globalSpaceIndex(spaceIdx, timeIdx);
return *materialParams_[globalIdx]; return *materialParams_[globalIdx];
} }
const MaterialLawParams &materialLawParams(unsigned globalIdx) const const MaterialLawParams& materialLawParams(unsigned globalIdx) const
{ return *materialParams_[globalIdx]; } { return *materialParams_[globalIdx]; }
/*! /*!
@@ -481,7 +481,9 @@ public:
* will need it one day? * will need it one day?
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return temperature_; } { return temperature_; }
// \} // \}
@@ -498,8 +500,10 @@ public:
* extraction and production wells, so all boundaries are no-flow. * extraction and production wells, so all boundaries are no-flow.
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values,
unsigned spaceIdx, unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
// no flow on top and bottom // no flow on top and bottom
values.setNoFlow(); values.setNoFlow();
@@ -519,7 +523,10 @@ public:
* the whole domain. * the whole domain.
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const void initial(PrimaryVariables& values,
const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
values.assignNaive(initialFluidState_); values.assignNaive(initialFluidState_);
@@ -538,13 +545,13 @@ public:
* saturated with a lower pressure than the remaining reservoir. * saturated with a lower pressure than the remaining reservoir.
*/ */
template <class Context> template <class Context>
void constraints(Constraints &constraints, const Context &context, void constraints(Constraints& constraints, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
if (this->simulator().episodeIndex() == 1) if (this->simulator().episodeIndex() == 1)
return; // no constraints during the "settle down" episode return; // no constraints during the "settle down" episode
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (isInjector_(pos)) { if (isInjector_(pos)) {
constraints.setActive(true); constraints.setActive(true);
constraints.assignNaive(injectorFluidState_); constraints.assignNaive(injectorFluidState_);
@@ -561,8 +568,10 @@ public:
* For this problem, the source term of all components is 0 everywhere. * For this problem, the source term of all components is 0 everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
@@ -570,7 +579,7 @@ public:
private: private:
void initFluidState_() void initFluidState_()
{ {
auto &fs = initialFluidState_; auto& fs = initialFluidState_;
////// //////
// set temperatures // set temperatures
@@ -590,7 +599,7 @@ private:
Scalar pw = pReservoir_; Scalar pw = pReservoir_;
PhaseVector pC; PhaseVector pC;
const auto &matParams = fineMaterialParams_; const auto& matParams = fineMaterialParams_;
MaterialLaw::capillaryPressures(pC, matParams, fs); MaterialLaw::capillaryPressures(pC, matParams, fs);
fs.setPressure(oilPhaseIdx, pw + (pC[oilPhaseIdx] - pC[waterPhaseIdx])); fs.setPressure(oilPhaseIdx, pw + (pC[oilPhaseIdx] - pC[waterPhaseIdx]));
@@ -643,8 +652,8 @@ private:
injFs.setSaturation(gasPhaseIdx, 0.0); injFs.setSaturation(gasPhaseIdx, 0.0);
// set the composition of the phases to immiscible // set the composition of the phases to immiscible
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
for (int compIdx = 0; compIdx < numComponents; ++compIdx) for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
injFs.setMoleFraction(phaseIdx, compIdx, 0.0); injFs.setMoleFraction(phaseIdx, compIdx, 0.0);
injFs.setMoleFraction(gasPhaseIdx, gasCompIdx, 1.0); injFs.setMoleFraction(gasPhaseIdx, gasCompIdx, 1.0);
@@ -676,7 +685,7 @@ private:
/*setEnthalpies=*/false); /*setEnthalpies=*/false);
} }
bool isProducer_(const GlobalPosition &pos) const bool isProducer_(const GlobalPosition& pos) const
{ {
Scalar x = pos[0] - this->boundingBoxMin()[0]; Scalar x = pos[0] - this->boundingBoxMin()[0];
Scalar y = pos[dim - 1] - this->boundingBoxMin()[dim - 1]; Scalar y = pos[dim - 1] - this->boundingBoxMin()[dim - 1];
@@ -691,7 +700,7 @@ private:
return width/2.0 - width*1e-5 < x && x < width/2.0 + width*(wellWidth_ + 1e-5); return width/2.0 - width*1e-5 < x && x < width/2.0 + width*(wellWidth_ + 1e-5);
} }
bool isInjector_(const GlobalPosition &pos) const bool isInjector_(const GlobalPosition& pos) const
{ {
Scalar x = pos[0] - this->boundingBoxMin()[0]; Scalar x = pos[0] - this->boundingBoxMin()[0];
Scalar y = pos[dim - 1] - this->boundingBoxMin()[dim - 1]; Scalar y = pos[dim - 1] - this->boundingBoxMin()[dim - 1];
@@ -706,7 +715,7 @@ private:
return x < width*wellWidth_ - width*1e-5 || x > width*(1.0 - wellWidth_) + width*1e-5; return x < width*wellWidth_ - width*1e-5 || x > width*(1.0 - wellWidth_) + width*1e-5;
} }
bool isFineMaterial_(const GlobalPosition &pos) const bool isFineMaterial_(const GlobalPosition& pos) const
{ return pos[dim - 1] > layerBottom_; } { return pos[dim - 1] > layerBottom_; }
DimMatrix fineK_; DimMatrix fineK_;

59
examples/problems/richardslensproblem.hh
View File

@@ -36,6 +36,7 @@
#include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp> #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
#include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp> #include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp>
#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp> #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -174,7 +175,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
RichardsLensProblem(Simulator &simulator) RichardsLensProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
, pnRef_(1e5) , pnRef_(1e5)
{ {
@@ -270,21 +271,21 @@ public:
* \copydoc FvBaseMultiPhaseProblem::temperature * \copydoc FvBaseMultiPhaseProblem::temperature
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ return temperature(context.globalSpaceIndex(spaceIdx, timeIdx), timeIdx); } { return temperature(context.globalSpaceIndex(spaceIdx, timeIdx), timeIdx); }
Scalar temperature(unsigned globalSpaceIdx, unsigned timeIdx) const Scalar temperature(unsigned OPM_UNUSED globalSpaceIdx, unsigned OPM_UNUSED timeIdx) const
{ return 273.15 + 10; } // -> 10°C { return 273.15 + 10; } // -> 10°C
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability * \copydoc FvBaseMultiPhaseProblem::intrinsicPermeability
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, const DimMatrix& intrinsicPermeability(const Context& context,
unsigned spaceIdx, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isInLens_(pos)) if (isInLens_(pos))
return lensK_; return lensK_;
return outerK_; return outerK_;
@@ -294,14 +295,16 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 0.4; } { return 0.4; }
/*! /*!
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& context,
unsigned spaceIdx, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
@@ -309,7 +312,8 @@ public:
return materialLawParams(globalSpaceIdx, timeIdx); return materialLawParams(globalSpaceIdx, timeIdx);
} }
const MaterialLawParams& materialLawParams(unsigned globalSpaceIdx, unsigned timeIdx) const const MaterialLawParams& materialLawParams(unsigned globalSpaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
if (dofIsInLens_[globalSpaceIdx]) if (dofIsInLens_[globalSpaceIdx])
return lensMaterialParams_; return lensMaterialParams_;
@@ -322,14 +326,15 @@ public:
* \copydetails Doxygen::contextParams * \copydetails Doxygen::contextParams
*/ */
template <class Context> template <class Context>
Scalar referencePressure(const Context &context, Scalar referencePressure(const Context& context,
unsigned spaceIdx, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ return referencePressure(context.globalSpaceIndex(spaceIdx, timeIdx), timeIdx); } { return referencePressure(context.globalSpaceIndex(spaceIdx, timeIdx), timeIdx); }
// the Richards model does not have an element context available at all places // the Richards model does not have an element context available at all places
// where the reference pressure is required... // where the reference pressure is required...
Scalar referencePressure(unsigned globalSpaceIdx, unsigned timeIdx) const Scalar referencePressure(unsigned OPM_UNUSED globalSpaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return pnRef_; } { return pnRef_; }
//! \} //! \}
@@ -343,15 +348,15 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, void boundary(BoundaryRateVector& values,
const Context &context, const Context& context,
unsigned spaceIdx, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(pos) || onRightBoundary_(pos)) { if (onLeftBoundary_(pos) || onRightBoundary_(pos)) {
const auto &materialParams = this->materialLawParams(context, spaceIdx, timeIdx); const auto& materialParams = this->materialLawParams(context, spaceIdx, timeIdx);
Scalar Sw = 0.0; Scalar Sw = 0.0;
Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs; Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs;
@@ -388,12 +393,12 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, void initial(PrimaryVariables& values,
const Context &context, const Context& context,
unsigned spaceIdx, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const auto &materialParams = this->materialLawParams(context, spaceIdx, timeIdx); const auto& materialParams = this->materialLawParams(context, spaceIdx, timeIdx);
Scalar Sw = 0.0; Scalar Sw = 0.0;
Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs; Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs;
@@ -412,35 +417,35 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, void source(RateVector& rate,
const Context &context, const Context& OPM_UNUSED context,
unsigned spaceIdx, unsigned OPM_UNUSED spaceIdx,
unsigned timeIdx) const unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < this->boundingBoxMin()[0] + eps_; } { return pos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < this->boundingBoxMin()[1] + eps_; } { return pos[1] < this->boundingBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool onInlet_(const GlobalPosition &pos) const bool onInlet_(const GlobalPosition& pos) const
{ {
Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0]; Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0];
Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width; Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width;
return onUpperBoundary_(pos) && 0.5 < lambda && lambda < 2.0 / 3.0; return onUpperBoundary_(pos) && 0.5 < lambda && lambda < 2.0 / 3.0;
} }
bool isInLens_(const GlobalPosition &pos) const bool isInLens_(const GlobalPosition& pos) const
{ {
for (unsigned i = 0; i < dimWorld; ++i) { for (unsigned i = 0; i < dimWorld; ++i) {
if (pos[i] < lensLowerLeft_[i] || pos[i] > lensUpperRight_[i]) if (pos[i] < lensLowerLeft_[i] || pos[i] > lensUpperRight_[i])

50
examples/problems/stokes2ctestproblem.hh
View File

@@ -28,7 +28,9 @@
#define EWOMS_STOKES_2C_TEST_PROBLEM_HH #define EWOMS_STOKES_2C_TEST_PROBLEM_HH
#include <ewoms/models/stokes/stokesmodel.hh> #include <ewoms/models/stokes/stokesmodel.hh>
#include <opm/material/fluidsystems/H2OAirFluidSystem.hpp> #include <opm/material/fluidsystems/H2OAirFluidSystem.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -126,7 +128,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
Stokes2cTestProblem(Simulator &simulator) Stokes2cTestProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -181,7 +183,9 @@ public:
* This problem assumes a temperature of 10 degrees Celsius. * This problem assumes a temperature of 10 degrees Celsius.
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 273.15 + 10; /* -> 10 deg C */ } { return 273.15 + 10; /* -> 10 deg C */ }
// \} // \}
@@ -199,21 +203,21 @@ public:
* upper edge. * upper edge.
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (onLowerBoundary_(pos)) if (onLowerBoundary_(pos))
values.setOutFlow(context, spaceIdx, timeIdx); values.setOutFlow(context, spaceIdx, timeIdx);
else if (onUpperBoundary_(pos)) { else if (onUpperBoundary_(pos))
// upper boundary is constraint! // upper boundary is constraint!
values = 0.0; values = 0.0;
} else
else {
// left and right boundaries // left and right boundaries
values.setNoFlow(context, spaceIdx, timeIdx); values.setNoFlow(context, spaceIdx, timeIdx);
}
} }
//! \} //! \}
@@ -231,10 +235,12 @@ public:
* 0.5% is set. * 0.5% is set.
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values,
const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &globalPos = context.pos(spaceIdx, timeIdx); const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
values = 0.0; values = 0.0;
// parabolic profile // parabolic profile
@@ -266,8 +272,10 @@ public:
* is 0 everywhere. * is 0 everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
/*! /*!
@@ -277,10 +285,12 @@ public:
* initial conditions. * initial conditions.
*/ */
template <class Context> template <class Context>
void constraints(Constraints &constraints, const Context &context, void constraints(Constraints& constraints,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onUpperBoundary_(pos)) { if (onUpperBoundary_(pos)) {
constraints.setActive(true); constraints.setActive(true);
@@ -290,16 +300,16 @@ public:
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &globalPos) const bool onLeftBoundary_(const GlobalPosition& globalPos) const
{ return globalPos[0] < this->boundingBoxMin()[0] + eps_; } { return globalPos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &globalPos) const bool onRightBoundary_(const GlobalPosition& globalPos) const
{ return globalPos[0] > this->boundingBoxMax()[0] - eps_; } { return globalPos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &globalPos) const bool onLowerBoundary_(const GlobalPosition& globalPos) const
{ return globalPos[1] < this->boundingBoxMin()[1] + eps_; } { return globalPos[1] < this->boundingBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &globalPos) const bool onUpperBoundary_(const GlobalPosition& globalPos) const
{ return globalPos[1] > this->boundingBoxMax()[1] - eps_; } { return globalPos[1] > this->boundingBoxMax()[1] - eps_; }
Scalar eps_; Scalar eps_;

38
examples/problems/stokesnitestproblem.hh
View File

@@ -29,7 +29,9 @@
#include <ewoms/models/stokes/stokesmodel.hh> #include <ewoms/models/stokes/stokesmodel.hh>
#include <ewoms/io/simplexgridmanager.hh> #include <ewoms/io/simplexgridmanager.hh>
#include <opm/material/fluidsystems/H2OAirFluidSystem.hpp> #include <opm/material/fluidsystems/H2OAirFluidSystem.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -138,7 +140,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
StokesNiTestProblem(Simulator &simulator) StokesNiTestProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -199,10 +201,10 @@ public:
* \copydoc FvBaseProblem::boundary * \copydoc FvBaseProblem::boundary
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (onUpperBoundary_(pos)) if (onUpperBoundary_(pos))
values.setOutFlow(context, spaceIdx, timeIdx); values.setOutFlow(context, spaceIdx, timeIdx);
@@ -227,10 +229,10 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
Scalar moleFrac[numComponents]; Scalar moleFrac[numComponents];
@@ -276,8 +278,10 @@ public:
* is 0 everywhere. * is 0 everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
/*! /*!
@@ -287,10 +291,12 @@ public:
* adjacent to the inlet. * adjacent to the inlet.
*/ */
template <class Context> template <class Context>
void constraints(Constraints &constraints, const Context &context, void constraints(Constraints& constraints,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onLowerBoundary_(pos) || onUpperBoundary_(pos)) { if (onLowerBoundary_(pos) || onUpperBoundary_(pos)) {
constraints.setActive(true); constraints.setActive(true);
@@ -301,25 +307,25 @@ public:
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < this->boundingBoxMin()[0] + eps_; } { return pos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < this->boundingBoxMin()[1] + eps_; } { return pos[1] < this->boundingBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool onBoundary_(const GlobalPosition &pos) const bool onBoundary_(const GlobalPosition& pos) const
{ {
return onLeftBoundary_(pos) || onRightBoundary_(pos) return onLeftBoundary_(pos) || onRightBoundary_(pos)
|| onLowerBoundary_(pos) || onUpperBoundary_(pos); || onLowerBoundary_(pos) || onUpperBoundary_(pos);
} }
bool inLens_(const GlobalPosition &pos) const bool inLens_(const GlobalPosition& pos) const
{ return pos[0] < 0.75 && pos[0] > 0.25 && pos[1] < 0.75 && pos[1] > 0.25; } { return pos[0] < 0.75 && pos[0] > 0.25 && pos[1] < 0.75 && pos[1] > 0.25; }
Scalar eps_; Scalar eps_;

39
examples/problems/stokestestproblem.hh
View File

@@ -32,6 +32,7 @@
#include <opm/material/fluidsystems/H2ON2FluidSystem.hpp> #include <opm/material/fluidsystems/H2ON2FluidSystem.hpp>
#include <opm/material/fluidsystems/GasPhase.hpp> #include <opm/material/fluidsystems/GasPhase.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -132,7 +133,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
StokesTestProblem(Simulator &simulator) StokesTestProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ eps_ = 1e-6; } { eps_ = 1e-6; }
@@ -175,7 +176,9 @@ public:
* This problem assumes a constant temperature of 10 degrees Celsius. * This problem assumes a constant temperature of 10 degrees Celsius.
*/ */
template <class Context> template <class Context>
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar temperature(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ return 273.15 + 10; } // -> 10 deg C { return 273.15 + 10; } // -> 10 deg C
//! \} //! \}
@@ -193,10 +196,10 @@ public:
* a parabolic velocity profile via constraints. * a parabolic velocity profile via constraints.
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
Scalar y = pos[1] - this->boundingBoxMin()[1]; Scalar y = pos[1] - this->boundingBoxMin()[1];
Scalar height = this->boundingBoxMax()[1] - this->boundingBoxMin()[1]; Scalar height = this->boundingBoxMax()[1] - this->boundingBoxMin()[1];
@@ -234,10 +237,10 @@ public:
* \copydoc FvBaseProblem::initial * \copydoc FvBaseProblem::initial
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
Scalar y = pos[1] - this->boundingBoxMin()[1]; Scalar y = pos[1] - this->boundingBoxMin()[1];
Scalar height = this->boundingBoxMax()[1] - this->boundingBoxMin()[1]; Scalar height = this->boundingBoxMax()[1] - this->boundingBoxMin()[1];
@@ -264,8 +267,10 @@ public:
* is 0 everywhere. * is 0 everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, const Context &context, unsigned spaceIdx, void source(RateVector& rate,
unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = Scalar(0.0); } { rate = Scalar(0.0); }
/*! /*!
@@ -275,10 +280,12 @@ public:
* velocity profile using constraints. * velocity profile using constraints.
*/ */
template <class Context> template <class Context>
void constraints(Constraints &constraints, const Context &context, void constraints(Constraints& constraints,
unsigned spaceIdx, unsigned timeIdx) const const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (onLeftBoundary_(pos) || onRightBoundary_(pos)) { if (onLeftBoundary_(pos) || onRightBoundary_(pos)) {
constraints.setActive(true); constraints.setActive(true);
@@ -289,19 +296,19 @@ public:
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < this->boundingBoxMin()[0] + eps_; } { return pos[0] < this->boundingBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < this->boundingBoxMin()[1] + eps_; } { return pos[1] < this->boundingBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool onBoundary_(const GlobalPosition &pos) const bool onBoundary_(const GlobalPosition& pos) const
{ {
return onLeftBoundary_(pos) || onRightBoundary_(pos) return onLeftBoundary_(pos) || onRightBoundary_(pos)
|| onLowerBoundary_(pos) || onUpperBoundary_(pos); || onLowerBoundary_(pos) || onUpperBoundary_(pos);

66
examples/problems/waterairproblem.hh
View File

@@ -39,6 +39,7 @@
#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp> #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
#include <opm/material/heatconduction/Somerton.hpp> #include <opm/material/heatconduction/Somerton.hpp>
#include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp> #include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/grid/yaspgrid.hh> #include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/dgfparser/dgfyasp.hh> #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
@@ -208,7 +209,7 @@ public:
/*! /*!
* \copydoc Doxygen::defaultProblemConstructor * \copydoc Doxygen::defaultProblemConstructor
*/ */
WaterAirProblem(Simulator &simulator) WaterAirProblem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
{ } { }
@@ -301,9 +302,9 @@ public:
* permeable than the lower one. * permeable than the lower one.
*/ */
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, unsigned spaceIdx, unsigned timeIdx) const const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineK_; return fineK_;
return coarseK_; return coarseK_;
@@ -313,9 +314,9 @@ public:
* \copydoc FvBaseMultiPhaseProblem::porosity * \copydoc FvBaseMultiPhaseProblem::porosity
*/ */
template <class Context> template <class Context>
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar porosity(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return finePorosity_; return finePorosity_;
else else
@@ -326,11 +327,11 @@ public:
* \copydoc FvBaseMultiPhaseProblem::materialLawParams * \copydoc FvBaseMultiPhaseProblem::materialLawParams
*/ */
template <class Context> template <class Context>
const MaterialLawParams& materialLawParams(const Context &context, const MaterialLawParams& materialLawParams(const Context& context,
unsigned spaceIdx, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineMaterialParams_; return fineMaterialParams_;
else else
@@ -343,7 +344,9 @@ public:
* In this case, we assume the rock-matrix to be granite. * In this case, we assume the rock-matrix to be granite.
*/ */
template <class Context> template <class Context>
Scalar heatCapacitySolid(const Context &context, unsigned spaceIdx, unsigned timeIdx) const Scalar heatCapacitySolid(const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ {
return return
790 // specific heat capacity of granite [J / (kg K)] 790 // specific heat capacity of granite [J / (kg K)]
@@ -355,9 +358,9 @@ public:
*/ */
template <class Context> template <class Context>
const HeatConductionLawParams& const HeatConductionLawParams&
heatConductionParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const heatConductionParams(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
if (isFineMaterial_(pos)) if (isFineMaterial_(pos))
return fineHeatCondParams_; return fineHeatCondParams_;
return coarseHeatCondParams_; return coarseHeatCondParams_;
@@ -379,11 +382,11 @@ public:
* right boundaries of the domain. * right boundaries of the domain.
*/ */
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, void boundary(BoundaryRateVector& values,
const Context &context, const Context& context,
unsigned spaceIdx, unsigned timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const auto &pos = context.cvCenter(spaceIdx, timeIdx); const auto& pos = context.cvCenter(spaceIdx, timeIdx);
assert(onLeftBoundary_(pos) || assert(onLeftBoundary_(pos) ||
onLowerBoundary_(pos) || onLowerBoundary_(pos) ||
onRightBoundary_(pos) || onRightBoundary_(pos) ||
@@ -432,12 +435,15 @@ public:
* liquid water and assume hydrostatic pressure. * liquid water and assume hydrostatic pressure.
*/ */
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const void initial(PrimaryVariables& values,
const Context& context,
unsigned spaceIdx,
unsigned timeIdx) const
{ {
Opm::CompositionalFluidState<Scalar, FluidSystem> fs; Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
initialFluidState_(fs, context, spaceIdx, timeIdx); initialFluidState_(fs, context, spaceIdx, timeIdx);
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
values.assignMassConservative(fs, matParams, /*inEquilibrium=*/true); values.assignMassConservative(fs, matParams, /*inEquilibrium=*/true);
} }
@@ -448,38 +454,40 @@ public:
* everywhere. * everywhere.
*/ */
template <class Context> template <class Context>
void source(RateVector &rate, void source(RateVector& rate,
const Context &context, unsigned spaceIdx, unsigned timeIdx) const const Context& OPM_UNUSED context,
unsigned OPM_UNUSED spaceIdx,
unsigned OPM_UNUSED timeIdx) const
{ rate = 0; } { rate = 0; }
//! \} //! \}
private: private:
bool onLeftBoundary_(const GlobalPosition &pos) const bool onLeftBoundary_(const GlobalPosition& pos) const
{ return pos[0] < eps_; } { return pos[0] < eps_; }
bool onRightBoundary_(const GlobalPosition &pos) const bool onRightBoundary_(const GlobalPosition& pos) const
{ return pos[0] > this->boundingBoxMax()[0] - eps_; } { return pos[0] > this->boundingBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &pos) const bool onLowerBoundary_(const GlobalPosition& pos) const
{ return pos[1] < eps_; } { return pos[1] < eps_; }
bool onUpperBoundary_(const GlobalPosition &pos) const bool onUpperBoundary_(const GlobalPosition& pos) const
{ return pos[1] > this->boundingBoxMax()[1] - eps_; } { return pos[1] > this->boundingBoxMax()[1] - eps_; }
bool onInlet_(const GlobalPosition &pos) const bool onInlet_(const GlobalPosition& pos) const
{ return onLowerBoundary_(pos) && (15.0 < pos[0]) && (pos[0] < 25.0); } { return onLowerBoundary_(pos) && (15.0 < pos[0]) && (pos[0] < 25.0); }
bool inHighTemperatureRegion_(const GlobalPosition &pos) const bool inHighTemperatureRegion_(const GlobalPosition& pos) const
{ return (20 < pos[0]) && (pos[0] < 30) && (pos[1] < 30); } { return (20 < pos[0]) && (pos[0] < 30) && (pos[1] < 30); }
template <class Context, class FluidState> template <class Context, class FluidState>
void initialFluidState_(FluidState &fs, void initialFluidState_(FluidState& fs,
const Context &context, const Context& context,
unsigned spaceIdx, unsigned spaceIdx,
unsigned timeIdx) const unsigned timeIdx) const
{ {
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx); const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
Scalar densityW = 1000.0; Scalar densityW = 1000.0;
fs.setPressure(liquidPhaseIdx, 1e5 + (maxDepth_ - pos[1])*densityW*9.81); fs.setPressure(liquidPhaseIdx, 1e5 + (maxDepth_ - pos[1])*densityW*9.81);
@@ -495,7 +503,7 @@ private:
// set the gas saturation and pressure // set the gas saturation and pressure
fs.setSaturation(gasPhaseIdx, 0); fs.setSaturation(gasPhaseIdx, 0);
Scalar pc[numPhases]; Scalar pc[numPhases];
const auto &matParams = materialLawParams(context, spaceIdx, timeIdx); const auto& matParams = materialLawParams(context, spaceIdx, timeIdx);
MaterialLaw::capillaryPressures(pc, matParams, fs); MaterialLaw::capillaryPressures(pc, matParams, fs);
fs.setPressure(gasPhaseIdx, fs.pressure(liquidPhaseIdx) + (pc[gasPhaseIdx] - pc[liquidPhaseIdx])); fs.setPressure(gasPhaseIdx, fs.pressure(liquidPhaseIdx) + (pc[gasPhaseIdx] - pc[liquidPhaseIdx]));
@@ -504,7 +512,7 @@ private:
CFRP::solve(fs, paramCache, liquidPhaseIdx, /*setViscosity=*/false, /*setEnthalpy=*/true); CFRP::solve(fs, paramCache, liquidPhaseIdx, /*setViscosity=*/false, /*setEnthalpy=*/true);
} }
void computeHeatCondParams_(HeatConductionLawParams &params, Scalar poro) void computeHeatCondParams_(HeatConductionLawParams& params, Scalar poro)
{ {
Scalar lambdaGranite = 2.8; // [W / (K m)] Scalar lambdaGranite = 2.8; // [W / (K m)]
@@ -538,7 +546,7 @@ private:
} }
} }
bool isFineMaterial_(const GlobalPosition &pos) const bool isFineMaterial_(const GlobalPosition& pos) const
{ return pos[dim-1] > layerBottom_; } { return pos[dim-1] > layerBottom_; }
DimMatrix fineK_; DimMatrix fineK_;

36
examples/tutorial1problem.hh
View File

@@ -167,7 +167,7 @@ class Tutorial1Problem
public: public:
//! The constructor of the problem. This only _allocates_ the memory required by the //! The constructor of the problem. This only _allocates_ the memory required by the
//! problem. The constructor is supposed to _never ever_ throw an exception. //! problem. The constructor is supposed to _never ever_ throw an exception.
Tutorial1Problem(Simulator &simulator) Tutorial1Problem(Simulator& simulator)
: ParentType(simulator) : ParentType(simulator)
, eps_(3e-6) , eps_(3e-6)
{ } { }
@@ -200,36 +200,37 @@ public:
//! Returns the temperature at a given position. //! Returns the temperature at a given position.
template <class Context> template <class Context>
Scalar temperature(const Context &context, int spaceIdx, int timeIdx) const Scalar temperature(const Context& /*context*/,
unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
{ return 283.15; } { return 283.15; }
//! Returns the intrinsic permeability tensor [m^2] at a position. //! Returns the intrinsic permeability tensor [m^2] at a position.
template <class Context> template <class Context>
const DimMatrix &intrinsicPermeability(const Context &context, /*@\label{tutorial1:permeability}@*/ const DimMatrix& intrinsicPermeability(const Context& /*context*/, /*@\label{tutorial1:permeability}@*/
int spaceIdx, int timeIdx) const unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
{ return K_; } { return K_; }
//! Defines the porosity [-] of the medium at a given position //! Defines the porosity [-] of the medium at a given position
template <class Context> template <class Context>
Scalar porosity(const Context &context, Scalar porosity(const Context& /*context*/,
int spaceIdx, int timeIdx) const /*@\label{tutorial1:porosity}@*/ unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const /*@\label{tutorial1:porosity}@*/
{ return 0.2; } { return 0.2; }
//! Returns the parameter object for the material law at a given position //! Returns the parameter object for the material law at a given position
template <class Context> template <class Context>
const MaterialLawParams &materialLawParams(const Context &context, /*@\label{tutorial1:matLawParams}@*/ const MaterialLawParams& materialLawParams(const Context& /*context*/, /*@\label{tutorial1:matLawParams}@*/
int spaceIdx, int timeIdx) const unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
{ return materialParams_; } { return materialParams_; }
//! Evaluates the boundary conditions. //! Evaluates the boundary conditions.
template <class Context> template <class Context>
void boundary(BoundaryRateVector &values, const Context &context, void boundary(BoundaryRateVector& values, const Context& context,
int spaceIdx, int timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
const auto &pos = context.pos(spaceIdx, timeIdx); const auto& pos = context.pos(spaceIdx, timeIdx);
if (pos[0] < eps_) { if (pos[0] < eps_) {
// Free-flow conditions on left boundary // Free-flow conditions on left boundary
const auto &materialParams = this->materialLawParams(context, spaceIdx, timeIdx); const auto& materialParams = this->materialLawParams(context, spaceIdx, timeIdx);
Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs; Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs;
Scalar Sw = 1.0; Scalar Sw = 1.0;
@@ -261,8 +262,8 @@ public:
//! position of the domain [kg/(m^3 * s)]. Positive values mean that //! position of the domain [kg/(m^3 * s)]. Positive values mean that
//! mass is created. //! mass is created.
template <class Context> template <class Context>
void source(RateVector &source, const Context &context, int spaceIdx, void source(RateVector& source, const Context& /*context*/,
int timeIdx) const unsigned /*spaceIdx*/, unsigned /*timeIdx*/) const
{ {
source[contiWettingEqIdx] = 0.0; source[contiWettingEqIdx] = 0.0;
source[contiNonWettingEqIdx] = 0.0; source[contiNonWettingEqIdx] = 0.0;
@@ -270,8 +271,8 @@ public:
//! Evaluates the initial value at a given position in the domain. //! Evaluates the initial value at a given position in the domain.
template <class Context> template <class Context>
void initial(PrimaryVariables &values, const Context &context, int spaceIdx, void initial(PrimaryVariables& values, const Context& context,
int timeIdx) const unsigned spaceIdx, unsigned timeIdx) const
{ {
Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs; Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs;
@@ -285,8 +286,7 @@ public:
// set pressure of the wetting phase to 200 kPa = 2 bar // set pressure of the wetting phase to 200 kPa = 2 bar
Scalar pC[numPhases]; Scalar pC[numPhases];
MaterialLaw::capillaryPressures(pC, materialLawParams(context, spaceIdx, MaterialLaw::capillaryPressures(pC, materialLawParams(context, spaceIdx, timeIdx),
timeIdx),
fs); fs);
fs.setPressure(wettingPhaseIdx, 200e3); fs.setPressure(wettingPhaseIdx, 200e3);
fs.setPressure(nonWettingPhaseIdx, 200e3 + pC[nonWettingPhaseIdx] - pC[nonWettingPhaseIdx]); fs.setPressure(nonWettingPhaseIdx, 200e3 + pC[nonWettingPhaseIdx] - pC[nonWettingPhaseIdx]);

20
tests/models/test_quadrature.cpp
View File

@@ -41,6 +41,14 @@
#include <ewoms/disc/vcfv/vcfvstencil.hh> #include <ewoms/disc/vcfv/vcfvstencil.hh>
#include <opm/material/common/Unused.hpp>
#if HAVE_DUNE_ALUGRID
#define EWOMS_NO_ALUGRID_UNUSED
#else
#define EWOMS_NO_ALUGRID_UNUSED OPM_UNUSED
#endif
const unsigned dim = 3; const unsigned dim = 3;
typedef double Scalar; typedef double Scalar;
typedef Ewoms::QuadrialteralQuadratureGeometry<Scalar, dim> QuadratureGeom; typedef Ewoms::QuadrialteralQuadratureGeometry<Scalar, dim> QuadratureGeom;
@@ -102,7 +110,7 @@ void testIdenityMapping()
} }
template <class Grid> template <class Grid>
void writeTetrahedronSubControlVolumes(const Grid &grid) void writeTetrahedronSubControlVolumes(const Grid& EWOMS_NO_ALUGRID_UNUSED grid)
{ {
#if HAVE_DUNE_ALUGRID #if HAVE_DUNE_ALUGRID
typedef typename Grid::LeafGridView GridView; typedef typename Grid::LeafGridView GridView;
@@ -134,14 +142,14 @@ void writeTetrahedronSubControlVolumes(const Grid &grid)
} }
} }
int cornerOffset = 0; unsigned cornerOffset = 0;
eIt = gridView.template begin<0>(); eIt = gridView.template begin<0>();
for (; eIt != eEndIt; ++eIt) { for (; eIt != eEndIt; ++eIt) {
stencil.update(*eIt); stencil.update(*eIt);
for (unsigned scvIdx = 0; scvIdx < stencil.numDof(); ++scvIdx) { for (unsigned scvIdx = 0; scvIdx < stencil.numDof(); ++scvIdx) {
const auto &scvLocalGeom = stencil.subControlVolume(scvIdx).localGeometry(); const auto &scvLocalGeom = stencil.subControlVolume(scvIdx).localGeometry();
std::vector<unsigned int> vertexIndices; std::vector<unsigned> vertexIndices;
for (unsigned i = 0; i < scvLocalGeom.numCorners; ++i) { for (unsigned i = 0; i < scvLocalGeom.numCorners; ++i) {
vertexIndices.push_back(cornerOffset); vertexIndices.push_back(cornerOffset);
++cornerOffset; ++cornerOffset;
@@ -185,7 +193,7 @@ void testTetrahedron()
} }
template <class Grid> template <class Grid>
void writeCubeSubControlVolumes(const Grid &grid) void writeCubeSubControlVolumes(const Grid& EWOMS_NO_ALUGRID_UNUSED grid)
{ {
#if HAVE_DUNE_ALUGRID #if HAVE_DUNE_ALUGRID
typedef typename Grid::LeafGridView GridView; typedef typename Grid::LeafGridView GridView;
@@ -215,7 +223,7 @@ void writeCubeSubControlVolumes(const Grid &grid)
} }
} }
int cornerOffset = 0; unsigned cornerOffset = 0;
eIt = gridView.template begin<0>(); eIt = gridView.template begin<0>();
for (; eIt != eEndIt; ++eIt) { for (; eIt != eEndIt; ++eIt) {
stencil.update(*eIt); stencil.update(*eIt);
@@ -318,7 +326,7 @@ void testQuadrature()
const auto &scvLocalGeom = stencil.subControlVolume(scvIdx).localGeometry(); const auto &scvLocalGeom = stencil.subControlVolume(scvIdx).localGeometry();
Dune::GeometryType geomType = scvLocalGeom.type(); Dune::GeometryType geomType = scvLocalGeom.type();
static const int quadratureOrder = 2; static const unsigned quadratureOrder = 2;
const auto &rule const auto &rule
= Dune::QuadratureRules<Scalar, dim>::rule(geomType, = Dune::QuadratureRules<Scalar, dim>::rule(geomType,
quadratureOrder); quadratureOrder);