diff --git a/appl/lecture/msm/buckleyleverett/buckleyleverett_analytic.hh b/appl/lecture/msm/buckleyleverett/buckleyleverett_analytic.hh
index 8f13316b4cf35fa1c68b296ab3642d83d5057bae..a65078a7c2a888644819ddbc8b39227f80fa64e4 100644
--- a/appl/lecture/msm/buckleyleverett/buckleyleverett_analytic.hh
+++ b/appl/lecture/msm/buckleyleverett/buckleyleverett_analytic.hh
@@ -111,8 +111,8 @@ public:
// std::cout<<"satVec_ = "<<satVec_<<std::endl;
FluidState fluidState;
- Scalar temp = problem_.temperature(dummyGlobal_, dummyElement_);
- Scalar press = problem_.referencePressure(dummyGlobal_, dummyElement_);
+ Scalar temp = problem_.temperature(dummyElement_);
+ Scalar press = problem_.referencePressure(dummyElement_);
Scalar viscosityW = FluidSystem::phaseViscosity(wPhaseIdx, temp, press, fluidState);
Scalar viscosityNW = FluidSystem::phaseViscosity(nPhaseIdx, temp, press, fluidState);
diff --git a/appl/lecture/msm/buckleyleverett/buckleyleverett_ff.cc b/appl/lecture/msm/buckleyleverett/buckleyleverett_ff.cc
index 95a2e349591035a5634189e32d40bd23558280a1..b8c9190eb3b741448187ba8c2a7107443a2fcfc4 100644
--- a/appl/lecture/msm/buckleyleverett/buckleyleverett_ff.cc
+++ b/appl/lecture/msm/buckleyleverett/buckleyleverett_ff.cc
@@ -59,6 +59,7 @@ int main(int argc, char** argv)
typedef GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
typedef GET_PROP_TYPE(TypeTag, PTAG(Problem)) Problem;
+ typedef GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
typedef Dune::FieldVector<Scalar, Grid::dimensionworld> GlobalPosition;
static const int dim = Grid::dimension;
@@ -109,10 +110,11 @@ int main(int argc, char** argv)
WettingPhase::Component::setDensity(interfaceFluidProps.IFP_DensityWettingFluid);
NonwettingPhase::Component::setDensity(interfaceFluidProps.IFP_DensityNonWettingFluid);
- Problem problem(grid.leafView(), lowerLeft, upperRight);
+ TimeManager timeManager;
+ Problem problem(timeManager, grid.leafView(), lowerLeft, upperRight);
- problem.timeManager().init(problem, 0, dt, tEnd);
- problem.timeManager().run();
+ timeManager.init(problem, 0, dt, tEnd);
+ timeManager.run();
return 0;
}
catch (Dune::Exception &e) {
diff --git a/appl/lecture/msm/buckleyleverett/buckleyleverettproblem.hh b/appl/lecture/msm/buckleyleverett/buckleyleverettproblem.hh
index 80ef62ab19a7e7cfc5ace139c08a2ae3c670d54d..c611d1c9c20bb7c0bf4cca58fa625f55a1e1a21d 100644
--- a/appl/lecture/msm/buckleyleverett/buckleyleverettproblem.hh
+++ b/appl/lecture/msm/buckleyleverett/buckleyleverettproblem.hh
@@ -125,10 +125,9 @@ SET_SCALAR_PROP(BuckleyLeverettProblem, CFLFactor, 0.95); //0.95
* \ingroup DecoupledProblems
*/
template<class TypeTag = TTAG(BuckleyLeverettProblem)>
-class BuckleyLeverettProblem: public IMPESProblem2P<TypeTag, BuckleyLeverettProblem<TypeTag> >
+class BuckleyLeverettProblem: public IMPESProblem2P<TypeTag>
{
- typedef BuckleyLeverettProblem<TypeTag> ThisType;
- typedef IMPESProblem2P<TypeTag, ThisType> ParentType;
+ typedef IMPESProblem2P<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
@@ -136,14 +135,22 @@ class BuckleyLeverettProblem: public IMPESProblem2P<TypeTag, BuckleyLeverettProb
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidState)) FluidState;
+ typedef typename GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
+
+ typedef typename GET_PROP_TYPE(TypeTag, PTAG(BoundaryTypes)) BoundaryTypes;
+ typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes))::PrimaryVariables PrimaryVariables;
+
enum
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
};
enum
{
- wetting = 0, nonwetting = 1
+ wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
+ eqIdxPress = Indices::pressureEq,
+ eqIdxSat = Indices::saturationEq
};
+
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef typename GridView::Traits::template Codim<0>::Entity Element;
@@ -151,11 +158,11 @@ class BuckleyLeverettProblem: public IMPESProblem2P<TypeTag, BuckleyLeverettProb
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
public:
- BuckleyLeverettProblem(const GridView &gridView,
+ BuckleyLeverettProblem(TimeManager& timeManager, const GridView &gridView,
const GlobalPosition lowerLeft = 0,
const GlobalPosition upperRight = 0,
const Scalar pleftbc = 2e5)
- : ParentType(gridView),
+ : ParentType(timeManager, gridView),
lowerLeft_(lowerLeft),
upperRight_(upperRight),
eps_(1e-8),
@@ -208,99 +215,78 @@ public:
analyticSolution_.addOutputVtkFields(this->resultWriter());
}
- //! Write the fields current solution into an VTK output file.
- void writeOutput()
- {
-// if (this->timeManager().time() > 43.1999e6)
-// {
- if (this->gridView().comm().rank() == 0)
- std::cout << "Writing result file for current time step\n";
-
- this->resultWriter().beginWrite(this->timeManager().time() + this->timeManager().timeStepSize());
- this->asImp_().addOutputVtkFields();
- this->resultWriter().endWrite();
-// }
- }
-
/*!
* \brief Returns the temperature within the domain.
*
* This problem assumes a temperature of 10 degrees Celsius.
*/
- Scalar temperature(const GlobalPosition& globalPos, const Element& element) const
+ Scalar temperatureAtPos(const GlobalPosition& globalPos) const
{
return 273.15 + 10; // -> 10°C
}
// \}
-
- Scalar referencePressure(const GlobalPosition& globalPos, const Element& element) const
+ Scalar referencePressureAtPos(const GlobalPosition& globalPos) const
{
- return 1e5;
- }
-
- std::vector<Scalar> source(const GlobalPosition& globalPos, const Element& element)
- {
- return std::vector<Scalar>(2, 0.0);
- }
-
- typename BoundaryConditions::Flags bctypePress(const GlobalPosition& globalPos, const Intersection& intersection) const
- {
- if (globalPos[0] < eps_)
- return BoundaryConditions::dirichlet;
-
- // all other boundaries
- else
- return BoundaryConditions::neumann;
+ return 1e5; // -> 10°C
}
- BoundaryConditions::Flags bctypeSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void sourceAtPos(PrimaryVariables &values,const GlobalPosition& globalPos) const
{
- if (globalPos[0] < eps_)
- return Dumux::BoundaryConditions::dirichlet;
- else if (globalPos[0] > upperRight_[0] - eps_)
- return Dumux::BoundaryConditions::outflow;
- else
- return Dumux::BoundaryConditions::neumann;
+ values = 0;
}
- Scalar dirichletPress(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos) const
{
- return pLeftBc_;
+ if (globalPos[0] < eps_)
+ {
+ bcTypes.setAllDirichlet();
+ }
+ else if (globalPos[0] > upperRight_[0] - eps_)
+ {
+ bcTypes.setNeumann(eqIdxPress);
+ bcTypes.setOutflow(eqIdxSat);
+ }
+ // all other boundaries
+ else
+ {
+ bcTypes.setAllNeumann();
+ }
}
- Scalar dirichletSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void dirichletAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
if (globalPos[0] < eps_)
- return 0.8;
-
- // all other boundaries
-
+ {
+ values[eqIdxPress] = pLeftBc_;
+ values[eqIdxSat] = 0.8;
+ }
else
- return 0.2;
+ {
+ values[eqIdxPress] = pLeftBc_;
+ values[eqIdxSat] = 0.2;
+ }
}
- std::vector<Scalar> neumann(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void neumannAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
- std::vector<Scalar> neumannFlux(2, 0.0);
+ values = 0;
if (globalPos[0]> upperRight_[0] - eps_)
{
// the volume flux should remain constant, when density is changed
// here, we multiply by the density of the NonWetting Phase
const Scalar referenceDensity = 1000.0;
- neumannFlux[nonwetting] = 3e-4 * densityNonWetting_/referenceDensity;
+
+ values[nPhaseIdx] = 3e-4 * densityNonWetting_/referenceDensity;
}
- return neumannFlux;
}
- Scalar initSat(const GlobalPosition& globalPos, const Element& element) const
+ void initialAtPos(PrimaryVariables &values,
+ const GlobalPosition &globalPos) const
{
- if (globalPos[0] < eps_)
- return 0.8;
-
- else
- return 0.2;
+ values[eqIdxPress] = 0;
+ values[eqIdxSat] = 0.2;
}
private:
diff --git a/appl/lecture/msm/mcwhorter/mcwhorter_analytic.hh b/appl/lecture/msm/mcwhorter/mcwhorter_analytic.hh
index c27c060c17f017241616fb02aef378649bbdbac6..4970225a061dac76b2a8be64e2f121929a21fa91 100644
--- a/appl/lecture/msm/mcwhorter/mcwhorter_analytic.hh
+++ b/appl/lecture/msm/mcwhorter/mcwhorter_analytic.hh
@@ -48,6 +48,8 @@ class McWhorterAnalytic
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidState)) FluidState;
+ typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes))::PrimaryVariables PrimaryVariables;
+
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
enum
@@ -56,7 +58,9 @@ class McWhorterAnalytic
};
enum
{
- wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
+ wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
+ eqIdxPress = Indices::pressureEq,
+ eqIdxSat = Indices::saturationEq
};
typedef Dune::BlockVector<Dune::FieldVector<Scalar, 1> > BlockVector;
@@ -128,7 +132,9 @@ public:
snr_ = materialLawParams.Snr();
porosity_ = problem_.spatialParameters().porosity(dummyGlobal_, dummyElement_);
permeability_ = problem_.spatialParameters().intrinsicPermeability(dummyGlobal_, dummyElement_)[0][0];
- sInit_ = problem_.initSat(dummyGlobal_, dummyElement_);
+ PrimaryVariables initVec;
+ problem_.initial(initVec, dummyElement_);
+ sInit_ = initVec[eqIdxSat];
Scalar s0 =(1 - snr_ - swr_);
time_=0;
@@ -142,8 +148,8 @@ public:
satVec_[i]=satVec_[i-1]+h_;
}
FluidState fluidState;
- Scalar temp = problem_.temperature(dummyGlobal_, dummyElement_);
- Scalar press = problem_.referencePressure(dummyGlobal_, dummyElement_);
+ Scalar temp = problem_.temperature(dummyElement_);
+ Scalar press = problem_.referencePressure(dummyElement_);
Scalar viscosityW = FluidSystem::phaseViscosity(wPhaseIdx, temp, press, fluidState);
Scalar viscosityNW = FluidSystem::phaseViscosity(nPhaseIdx, temp, press, fluidState);
diff --git a/appl/lecture/msm/mcwhorter/mcwhorter_ff.cc b/appl/lecture/msm/mcwhorter/mcwhorter_ff.cc
index 62bbabf4de1e87364eb9f73cd7ef36a151b89350..b3a6068c43beeb36f884562c5ad11cd16cfb813b 100644
--- a/appl/lecture/msm/mcwhorter/mcwhorter_ff.cc
+++ b/appl/lecture/msm/mcwhorter/mcwhorter_ff.cc
@@ -30,6 +30,7 @@ int main(int argc, char** argv)
typedef GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
typedef GET_PROP_TYPE(TypeTag, PTAG(Problem)) Problem;
+ typedef GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
typedef Dune::FieldVector<Scalar, Grid::dimensionworld> GlobalPosition;
static const int dim = Grid::dimension;
@@ -78,10 +79,11 @@ int main(int argc, char** argv)
WettingPhase::Component::setViscosity(interfaceFluidProps.IFP_ViscosityWettingFluid);
NonwettingPhase::Component::setViscosity(interfaceFluidProps.IFP_ViscosityNonWettingFluid);
- Problem problem(grid.leafView(), lowerLeft, upperRight);
+ TimeManager timeManager;
+ Problem problem(timeManager, grid.leafView(), lowerLeft, upperRight);
- problem.timeManager().init(problem, 0, dt, tEnd);
- problem.timeManager().run();
+ timeManager.init(problem, 0, dt, tEnd);
+ timeManager.run();
return 0;
diff --git a/appl/lecture/msm/mcwhorter/mcwhorterproblem.hh b/appl/lecture/msm/mcwhorter/mcwhorterproblem.hh
index e07bd42b428c5aefa29286f399b6315316c6c339..4e70b431b784d939c3729828377d35385574d738 100644
--- a/appl/lecture/msm/mcwhorter/mcwhorterproblem.hh
+++ b/appl/lecture/msm/mcwhorter/mcwhorterproblem.hh
@@ -119,10 +119,9 @@ SET_SCALAR_PROP(McWhorterProblem, CFLFactor, 0.8);
//! @brief McWhorter transport problem
template<class TypeTag = TTAG(McWhorterProblem)>
-class McWhorterProblem: public IMPESProblem2P<TypeTag, McWhorterProblem<TypeTag> >
+class McWhorterProblem: public IMPESProblem2P<TypeTag>
{
- typedef McWhorterProblem<TypeTag> ThisType;
- typedef IMPESProblem2P<TypeTag, ThisType> ParentType;
+ typedef IMPESProblem2P<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
@@ -130,14 +129,22 @@ class McWhorterProblem: public IMPESProblem2P<TypeTag, McWhorterProblem<TypeTag>
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidState)) FluidState;
+ typedef typename GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
+
+ typedef typename GET_PROP_TYPE(TypeTag, PTAG(BoundaryTypes)) BoundaryTypes;
+ typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes))::PrimaryVariables PrimaryVariables;
+
enum
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
};
enum
{
- wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
+ wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
+ eqIdxPress = Indices::pressureEq,
+ eqIdxSat = Indices::saturationEq
};
+
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef typename GridView::Traits::template Codim<0>::Entity Element;
@@ -147,17 +154,19 @@ class McWhorterProblem: public IMPESProblem2P<TypeTag, McWhorterProblem<TypeTag>
public:
- McWhorterProblem(const GridView &gridView,
+ McWhorterProblem(TimeManager& timeManager, const GridView &gridView,
const GlobalPosition lowerLeft = 0,
const GlobalPosition upperRight = 0,
const Scalar pleftbc = 2.0e5)//? initial oil pressure?
- : ParentType(gridView),
+ : ParentType(timeManager, gridView),
lowerLeft_(lowerLeft),
upperRight_(upperRight),
eps_(1e-6),
pLeftBc_(pleftbc),
analyticSolution_(*this)
- {}
+ {
+ this->setOutputInterval(10);
+ }
/*!
* \brief The problem name.
@@ -187,73 +196,62 @@ public:
analyticSolution_.addOutputVtkFields(this->resultWriter());
}
- bool shouldWriteOutput() const
- {
- if (this->timeManager().timeStepIndex() % 10 == 0)
- {
- return true;
- }
- return false;
- }
-
/*!
* \brief Returns the temperature within the domain.
*
* This problem assumes a temperature of 10 degrees Celsius.
*/
-
- Scalar temperature(const GlobalPosition& globalPos, const Element& element) const
- {return 273.15 + 10; // -> 10°C
- }
-
- Scalar referencePressure(const GlobalPosition& globalPos, const Element& element) const
- {
- return 1e5;
- }
- std::vector<Scalar> source (const GlobalPosition& globalPos, const Element& element)
+ Scalar temperatureAtPos(const GlobalPosition& globalPos) const
{
- return std::vector<Scalar>(2,0.0);//no source and sink terms
+ return 273.15 + 10; // -> 10°C
}
- BoundaryConditions::Flags bctypePress(const GlobalPosition& globalPos, const Intersection& intersection) const
+ Scalar referencePressureAtPos(const GlobalPosition& globalPos) const
{
- if (globalPos[0] < eps_)//west
- return BoundaryConditions::dirichlet;
-
- // all other boundaries
- else
- return BoundaryConditions::neumann;
+ return 1e5; // -> 10°C
}
- BoundaryConditions::Flags bctypeSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void sourceAtPos(PrimaryVariables &values,const GlobalPosition& globalPos) const
{
- if (globalPos[0] < eps_)// west and east
- return Dumux::BoundaryConditions::dirichlet;
- else
- return Dumux::BoundaryConditions::neumann;
+ values = 0;
}
- Scalar dirichletPress(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos) const
{
- return pLeftBc_;
+ if (globalPos[0] < eps_)//west
+ {
+ bcTypes.setAllDirichlet();
+ }
+ // all other boundaries
+ else
+ {
+ bcTypes.setAllNeumann();
+ }
}
- Scalar dirichletSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void dirichletAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
- if (globalPos[0] < eps_)//west
- return 1.0;
- else //east, north and south are Neumann
- return 0.0;
+ if (globalPos[0] < eps_)
+ {
+ values[eqIdxPress] = pLeftBc_;
+ values[eqIdxSat] = 1.0;
+ }
+ else
+ {
+ values[eqIdxPress] = pLeftBc_;
+ values[eqIdxSat] = 0.0;
+ }
}
- std::vector<Scalar> neumann(const GlobalPosition& globalPos, const Intersection& intersection) const
+ void neumannAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
- return std::vector<Scalar>(2,0.0);
+ values = 0;
}
- Scalar initSat (const GlobalPosition& globalPos, const Element& element) const
+ void initialAtPos(PrimaryVariables &values,
+ const GlobalPosition &globalPos) const
{
- return 0.0;
+ values = 0;
}
/* McWhorterProblem(VC& variables, Fluid& wettingphase, Fluid& nonwettingphase,
diff --git a/test/decoupled/1p/benchmarkresult.hh b/test/decoupled/1p/benchmarkresult.hh
index 5c2ae389f4a1ce48dcbf60603a904fd739fa94ad..d09a89ea21a9f6edcee41a181b97b7f28f63031d 100644
--- a/test/decoupled/1p/benchmarkresult.hh
+++ b/test/decoupled/1p/benchmarkresult.hh
@@ -431,7 +431,9 @@ public:
uMin = std::min(uMin, approxPressure);
uMax = std::max(uMax, approxPressure);
- sumf += volume*(problem.source(global, element)[0]);
+ typename Problem::PrimaryVariables sourceVec;
+ problem.source(sourceVec, element);
+ sumf += volume*sourceVec[0];
// get the absolute permeability
Dune::FieldMatrix<double,dim,dim> K = problem.spatialParameters().intrinsicPermeability(global, element);
diff --git a/test/decoupled/1p/test_diffusion.cc b/test/decoupled/1p/test_diffusion.cc
index 1fd8f1daa3897ae782ce03cd42df4e508e3ffa7c..f6689a34aff4abfa72d5f9339613ac5973eab1b7 100644
--- a/test/decoupled/1p/test_diffusion.cc
+++ b/test/decoupled/1p/test_diffusion.cc
@@ -87,6 +87,7 @@ int main(int argc, char** argv)
typedef GET_PROP_TYPE(TTAG(FVVelocity2PTestProblem), PTAG(Problem)) FVProblem;
FVProblem fvProblem(grid.leafView(), delta);
+ fvProblem.setName("fvdiffusion");
timer.reset();
fvProblem.init();
fvProblem.model().calculateVelocity();
@@ -97,6 +98,7 @@ int main(int argc, char** argv)
typedef GET_PROP_TYPE(TTAG(FVMPFAOVelocity2PTestProblem), PTAG(Problem)) MPFAOProblem;
MPFAOProblem mpfaProblem(grid.leafView(), delta);
+ mpfaProblem.setName("fvmpfaodiffusion");
timer.reset();
mpfaProblem.init();
mpfaProblem.model().calculateVelocity();
@@ -107,6 +109,7 @@ int main(int argc, char** argv)
typedef GET_PROP_TYPE(TTAG(MimeticPressure2PTestProblem), PTAG(Problem)) MimeticProblem;
MimeticProblem mimeticProblem(grid.leafView(), delta);
+ mimeticProblem.setName("mimeticdiffusion");
timer.reset();
mimeticProblem.init();
mimeticProblem.model().calculateVelocity();
diff --git a/test/decoupled/1p/test_diffusion_problem.hh b/test/decoupled/1p/test_diffusion_problem.hh
index 22a2fdf84fa837fecf6a9b6472b24a750b581849..eeb08d154678dfd2eec4d346383575f263fb2a78 100644
--- a/test/decoupled/1p/test_diffusion_problem.hh
+++ b/test/decoupled/1p/test_diffusion_problem.hh
@@ -46,16 +46,6 @@
namespace Dumux
{
-//! \cond INTERNAL
-struct FileNames
-{
- enum
- {
- TPFAName = 0, MPFAName = 1, MimeticName = 2
- };
-};
-//! \endcond
-
/*!
* \ingroup IMPETtests
*/
@@ -123,25 +113,20 @@ public:
// Enable gravity
SET_BOOL_PROP(DiffusionTestProblem, EnableGravity, false);
-NEW_PROP_TAG( FileName );
-
// set the types for the 2PFA FV method
NEW_TYPE_TAG(FVVelocity2PTestProblem, INHERITS_FROM(DiffusionTestProblem));
-SET_INT_PROP(FVVelocity2PTestProblem, FileName, FileNames::TPFAName);
SET_INT_PROP(FVVelocity2PTestProblem, IterationNumberPreconditioner, 0);
SET_TYPE_PROP(FVVelocity2PTestProblem, Model, Dumux::FVVelocity2P<TTAG(FVVelocity2PTestProblem)>);
SET_TYPE_PROP(FVVelocity2PTestProblem, Problem, Dumux::TestDiffusionProblem<TTAG(FVVelocity2PTestProblem)>);
// set the types for the MPFA-O FV method
NEW_TYPE_TAG(FVMPFAOVelocity2PTestProblem, INHERITS_FROM(DiffusionTestProblem));
-SET_INT_PROP(FVMPFAOVelocity2PTestProblem, FileName, FileNames::MPFAName);
SET_INT_PROP(FVMPFAOVelocity2PTestProblem, IterationNumberPreconditioner, 1);
SET_TYPE_PROP(FVMPFAOVelocity2PTestProblem, Model, Dumux::FVMPFAOVelocity2P<TypeTag>);
SET_TYPE_PROP(FVMPFAOVelocity2PTestProblem, Problem, Dumux::TestDiffusionProblem<TTAG(FVMPFAOVelocity2PTestProblem)>);
// set the types for the mimetic FD method
-NEW_TYPE_TAG(MimeticPressure2PTestProblem, INHERITS_FROM(DiffusionTestProblem));
-SET_INT_PROP(MimeticPressure2PTestProblem, FileName, FileNames::MimeticName);
+NEW_TYPE_TAG(MimeticPressure2PTestProblem, INHERITS_FROM(DiffusionTestProblem, Mimetic));
SET_INT_PROP(MimeticPressure2PTestProblem, IterationNumberPreconditioner, 0);
SET_TYPE_PROP(MimeticPressure2PTestProblem, Model, Dumux::MimeticPressure2P<TTAG(MimeticPressure2PTestProblem)>);
SET_TYPE_PROP(MimeticPressure2PTestProblem, Problem, Dumux::TestDiffusionProblem<TTAG(MimeticPressure2PTestProblem)>);
@@ -153,10 +138,9 @@ SET_TYPE_PROP(MimeticPressure2PTestProblem, Problem, Dumux::TestDiffusionProblem
* \brief test problem for diffusion models from the FVCA5 benchmark.
*/
template<class TypeTag = TTAG(DiffusionTestProblem)>
-class TestDiffusionProblem: public DiffusionProblem2P<TypeTag, TestDiffusionProblem<TypeTag> >
+class TestDiffusionProblem: public DiffusionProblem2P<TypeTag>
{
- typedef TestDiffusionProblem<TypeTag> ThisType;
- typedef DiffusionProblem2P<TypeTag, ThisType> ParentType;
+ typedef DiffusionProblem2P<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
@@ -164,6 +148,8 @@ class TestDiffusionProblem: public DiffusionProblem2P<TypeTag, TestDiffusionProb
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidState)) FluidState;
+ typedef typename GET_PROP_TYPE(TypeTag, PTAG(BoundaryTypes)) BoundaryTypes;
+
enum
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
@@ -171,7 +157,9 @@ class TestDiffusionProblem: public DiffusionProblem2P<TypeTag, TestDiffusionProb
enum
{
- wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
+ wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
+ eqIdxPress = Indices::pressureEq,
+ eqIdxSat = Indices::saturationEq
};
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
@@ -179,9 +167,10 @@ class TestDiffusionProblem: public DiffusionProblem2P<TypeTag, TestDiffusionProb
typedef typename GridView::Traits::template Codim<0>::Entity Element;
typedef typename GridView::Intersection Intersection;
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
- typedef Dune::FieldVector<Scalar, dim> LocalPosition;
public:
+ typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes))::PrimaryVariables PrimaryVariables;
+
TestDiffusionProblem(const GridView &gridView, const double delta = 1.0) :
ParentType(gridView), delta_(delta)
{
@@ -194,47 +183,30 @@ public:
*/
// \{
- /*!
- * \brief The problem name.
- *
- * This is used as a prefix for files generated by the simulation.
- */
- const char *name() const
- {
- switch (GET_PROP_VALUE(TypeTag, PTAG(FileName)))
- {
- case FileNames::TPFAName:
- return "fvdiffusion";
- case FileNames::MPFAName:
- return "fvmpfaodiffusion";
- case FileNames::MimeticName:
- return "mimeticdiffusion";
- }
- return "test_diffusion";
- }
-
bool shouldWriteRestartFile() const
{ return false; }
/*!
- * \brief Returns the temperature within the domain.
- *
- * This problem assumes a temperature of 10 degrees Celsius.
- */
- Scalar temperature(const GlobalPosition& globalPos, const Element& element) const
+ * \brief Returns the temperature within the domain.
+ *
+ * This problem assumes a temperature of 10 degrees Celsius.
+ */
+ Scalar temperatureAtPos(const GlobalPosition& globalPos) const
{
return 273.15 + 10; // -> 10°C
}
// \}
- Scalar referencePressure(const GlobalPosition& globalPos, const Element& element) const
+
+ //! Returns the reference pressure for evaluation of constitutive relations
+ Scalar referencePressureAtPos(const GlobalPosition& globalPos) const
{
return 1e5; // -> 10°C
}
- std::vector<Scalar> source(const GlobalPosition& globalPos, const Element& element)
- {
+ void sourceAtPos(PrimaryVariables &values,const GlobalPosition& globalPos) const
+ {
double pi = 4.0*atan(1.0);
double rt = globalPos[0]*globalPos[0]+globalPos[1]*globalPos[1];
double ux = pi*cos(pi*globalPos[0])*sin(pi*globalPos[1]);
@@ -244,42 +216,36 @@ public:
double kyy = (globalPos[0]*globalPos[0] + delta_*globalPos[1]*globalPos[1])/rt;
double f0 = sin(pi*globalPos[0])*sin(pi*globalPos[1])*pi*pi*(1.0 + delta_)*(globalPos[0]*globalPos[0] + globalPos[1]*globalPos[1])
+ cos(pi*globalPos[0])*sin(pi*globalPos[1])*pi*(1.0 - 3.0*delta_)*globalPos[0]
- + cos(pi*globalPos[1])*sin(pi*globalPos[0])*pi*(1.0 - 3.0*delta_)*globalPos[1]
- + cos(pi*globalPos[1])*cos(pi*globalPos[0])*2.0*pi*pi*(1.0 - delta_)*globalPos[0]*globalPos[1];
-
- std::vector<double> result(2, 0.0);
- result[wPhaseIdx]=(f0 + 2.0*(globalPos[0]*(kxx*ux + kxy*uy) + globalPos[1]*(kxy*ux + kyy*uy)))/rt;
+ + cos(pi*globalPos[1])*sin(pi*globalPos[0])*pi*(1.0 - 3.0*delta_)*globalPos[1] + cos(pi*globalPos[1])*cos(pi*globalPos[0])*2.0*pi*pi*(1.0 - delta_)*globalPos[0]*globalPos[1];
- return (result);
- }
-
- typename BoundaryConditions::Flags bctypePress(const GlobalPosition& globalPos, const Intersection& intersection) const
- {
- return BoundaryConditions::dirichlet;
+ values[wPhaseIdx] = (f0 + 2.0*(globalPos[0]*(kxx*ux + kxy*uy) + globalPos[1]*(kxy*ux + kyy*uy)))/rt;
}
- Scalar dirichletPress(const GlobalPosition& globalPos, const Intersection& intersection) const
+ /*!
+ * \brief Returns the type of boundary condition.
+ *
+ *
+ * BC for saturation equation can be dirichlet (saturation), neumann (flux), or outflow.
+ */
+ void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos) const
{
- return (exact(globalPos));
+
+ bcTypes.setAllDirichlet();
}
- typename BoundaryConditions::Flags bctypeSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ //! set dirichlet condition (saturation [-])
+ void dirichletAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
- return BoundaryConditions::dirichlet;
+ values[eqIdxPress] = exact(globalPos);
+ values[eqIdxSat] = 1.0;
}
- Scalar dirichletSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ //! set neumann condition for phases (flux, [kg/(m^2 s)])
+ void neumannAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
- return 1.0;
+ values = 0;
}
- std::vector<Scalar> neumann(const GlobalPosition& globalPos, const Intersection& intersection) const
- {
- std::vector<Scalar> neumannFlux(2, 0.0);
-
- return neumannFlux;
- }
-
Scalar exact (const GlobalPosition& globalPos) const
{
double pi = 4.0*atan(1.0);
diff --git a/test/decoupled/2p/test_impes.cc b/test/decoupled/2p/test_impes.cc
index 15faa88e818b3889d62bef2b38b6a445cdd8fa97..e1acad1a056b80b224a8217134ccc29feb7479f6 100644
--- a/test/decoupled/2p/test_impes.cc
+++ b/test/decoupled/2p/test_impes.cc
@@ -54,6 +54,7 @@ int main(int argc, char** argv)
typedef GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
typedef GET_PROP_TYPE(TypeTag, PTAG(Problem)) Problem;
+ typedef GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
typedef Dune::FieldVector<Scalar, Grid::dimensionworld> GlobalPosition;
static const int dim = Grid::dimension;
@@ -106,14 +107,15 @@ int main(int argc, char** argv)
////////////////////////////////////////////////////////////
// instantiate and run the concrete problem
////////////////////////////////////////////////////////////
- Problem problem(grid.leafView(), lowerLeft, upperRight);
+ TimeManager timeManager;
+ Problem problem(timeManager, grid.leafView());
// load restart file if necessarry
if (restart)
problem.deserialize(restartTime);
- problem.timeManager().init(problem, 0, dt, tEnd, !restart);
- problem.timeManager().run();
+ timeManager.init(problem, 0, dt, tEnd, !restart);
+ timeManager.run();
return 0;
}
catch (Dune::Exception &e) {
diff --git a/test/decoupled/2p/test_impes_problem.hh b/test/decoupled/2p/test_impes_problem.hh
index a22fbb687662d1e7c13e0b5d154aa8abd7b3f11d..245484a4aa199cdbf3c16ec06c7fee793c236ca4 100644
--- a/test/decoupled/2p/test_impes_problem.hh
+++ b/test/decoupled/2p/test_impes_problem.hh
@@ -149,10 +149,9 @@ SET_SCALAR_PROP(IMPESTestProblem, CFLFactor, 0.95);
* where the argument defines the simulation endtime.
*/
template<class TypeTag = TTAG(IMPESTestProblem)>
-class TestIMPESProblem: public IMPESProblem2P<TypeTag, TestIMPESProblem<TypeTag> >
+class TestIMPESProblem: public IMPESProblem2P<TypeTag>
{
-typedef TestIMPESProblem<TypeTag> ThisType;
-typedef IMPESProblem2P<TypeTag, ThisType> ParentType;
+typedef IMPESProblem2P<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
@@ -160,6 +159,8 @@ typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidState)) FluidState;
+typedef typename GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
+
enum
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
@@ -167,7 +168,9 @@ enum
enum
{
- wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
+ wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
+ eqIdxPress = Indices::pressureEq,
+ eqIdxSat = Indices::saturationEq
};
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
@@ -175,11 +178,13 @@ typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef typename GridView::Traits::template Codim<0>::Entity Element;
typedef typename GridView::Intersection Intersection;
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
-typedef Dune::FieldVector<Scalar, dim> LocalPosition;
+
+typedef typename GET_PROP_TYPE(TypeTag, PTAG(BoundaryTypes)) BoundaryTypes;
+typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes))::PrimaryVariables PrimaryVariables;
public:
-TestIMPESProblem(const GridView &gridView, const GlobalPosition lowerLeft = 0, const GlobalPosition upperRight = 0) :
-ParentType(gridView), lowerLeft_(lowerLeft), upperRight_(upperRight)
+TestIMPESProblem(TimeManager& timeManager, const GridView &gridView) :
+ParentType(timeManager, gridView)
{
}
@@ -208,7 +213,7 @@ bool shouldWriteRestartFile() const
*
* This problem assumes a temperature of 10 degrees Celsius.
*/
-Scalar temperature(const GlobalPosition& globalPos, const Element& element) const
+Scalar temperatureAtPos(const GlobalPosition& globalPos) const
{
return 273.15 + 10; // -> 10°C
}
@@ -216,98 +221,88 @@ Scalar temperature(const GlobalPosition& globalPos, const Element& element) cons
// \}
//! Returns the reference pressure for evaluation of constitutive relations
-Scalar referencePressure(const GlobalPosition& globalPos, const Element& element) const
+Scalar referencePressureAtPos(const GlobalPosition& globalPos) const
{
return 1e5; // -> 10°C
}
-//!source term [kg/(m^3 s)]
-std::vector<Scalar> source(const GlobalPosition& globalPos, const Element& element)
+void sourceAtPos(PrimaryVariables &values,const GlobalPosition& globalPos) const
{
- return std::vector<Scalar>(2, 0.0);
+ values = 0;
}
/*!
-* \brief Returns the type of boundary condition for the pressure equation.
+* \brief Returns the type of boundary condition.
*
-* BC can be dirichlet (pressure) or neumann (flux).
-*/
-typename BoundaryConditions::Flags bctypePress(const GlobalPosition& globalPos, const Intersection& intersection) const
-{
- if ((globalPos[0] < eps_))
- return BoundaryConditions::dirichlet;
- // all other boundaries
- return BoundaryConditions::neumann;
-}
-
-/*!
-* \brief Returns the type of boundary condition for the saturation equation.
+* BC for pressure equation can be dirichlet (pressure) or neumann (flux).
*
-* BC can be dirichlet (saturation), neumann (flux), or outflow.
+* BC for saturation equation can be dirichlet (saturation), neumann (flux), or outflow.
*/
-BoundaryConditions::Flags bctypeSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos) const
{
- if (globalPos[0] < eps_)
- return Dumux::BoundaryConditions::dirichlet;
- else if (globalPos[0] > upperRight_[0] - eps_)
- return Dumux::BoundaryConditions::outflow;
- else
- return Dumux::BoundaryConditions::neumann;
+ if (globalPos[0] < eps_)
+ {
+ bcTypes.setAllDirichlet();
+ }
+ else if (globalPos[0] > this->bboxMax()[0] - eps_)
+ {
+ bcTypes.setNeumann(eqIdxPress);
+ bcTypes.setOutflow(eqIdxSat);
+ }
+ // all other boundaries
+ else
+ {
+ bcTypes.setAllNeumann();
+ }
}
-//! return dirichlet condition (pressure, [Pa])
-Scalar dirichletPress(const GlobalPosition& globalPos, const Intersection& intersection) const
+//! set dirichlet condition (pressure [Pa], saturation [-])
+void dirichletAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
+ values = 0;
if (globalPos[0] < eps_)
{
if (GET_PROP_VALUE(TypeTag, PTAG(EnableGravity)))
{
- const Element& element = *(intersection.inside());
-
- Scalar pRef = referencePressure(globalPos, element);
- Scalar temp = temperature(globalPos, element);
- Scalar sat = this->variables().satElement(element);
+ Scalar pRef = referencePressureAtPos(globalPos);
+ Scalar temp = temperatureAtPos(globalPos);
+ Scalar sat = 1;
FluidState fluidState;
fluidState.update(sat, pRef, pRef, temp);
- return (2e5 + (upperRight_[dim-1] - globalPos[dim-1]) * FluidSystem::phaseDensity(wPhaseIdx, temp, pRef, fluidState) * this->gravity().two_norm());
+ values[eqIdxPress] = (2e5 + (this->bboxMax()[dim-1] - globalPos[dim-1]) * FluidSystem::phaseDensity(wPhaseIdx, temp, pRef, fluidState) * this->gravity().two_norm());
}
else
- return 2e5;
+ {
+ values[eqIdxPress] = 2e5;
+ }
+ values[eqIdxSat] = 0.8;
+ }
+ else
+ {
+ values[eqIdxPress] = 2e5;
+ values[eqIdxSat] = 0.2;
}
- // all other boundaries
- return 2e5;
-}
-
-//! return dirichlet condition (saturation, [-])
-Scalar dirichletSat(const GlobalPosition& globalPos, const Intersection& intersection) const
-{
- if (globalPos[0] < eps_)
- return 0.8;
- // all other boundaries
- return 0.2;
}
-//! return neumann condition (flux, [kg/(m^2 s)])
-std::vector<Scalar> neumann(const GlobalPosition& globalPos, const Intersection& intersection) const
+//! set neumann condition for phases (flux, [kg/(m^2 s)])
+void neumannAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
- std::vector<Scalar> neumannFlux(2, 0.0);
- if (globalPos[0] > upperRight_[0] - eps_)
+ values = 0;
+ if (globalPos[0] > this->bboxMax()[0] - eps_)
{
- neumannFlux[nPhaseIdx] = 3e-4;
+ values[nPhaseIdx] = 3e-4;
}
- return neumannFlux;
}
-
-//! initial condition for saturation
-Scalar initSat(const GlobalPosition& globalPos, const Element& element) const
+//! return initial solution -> only saturation values have to be given!
+void initialAtPos(PrimaryVariables &values,
+ const GlobalPosition &globalPos) const
{
- return 0.2;
+ values[eqIdxPress] = 0;
+ values[eqIdxSat] = 0.2;
}
private:
-GlobalPosition lowerLeft_;
-GlobalPosition upperRight_;
static constexpr Scalar eps_ = 1e-6;
};
diff --git a/test/decoupled/2p/test_transport.cc b/test/decoupled/2p/test_transport.cc
index 8979ef98f3afa85b0512ed10e35c5adcc664ea37..c7132a053357ed9ccb9c014177b3f0a7a721bcd0 100644
--- a/test/decoupled/2p/test_transport.cc
+++ b/test/decoupled/2p/test_transport.cc
@@ -54,6 +54,7 @@ int main(int argc, char** argv)
typedef GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
typedef GET_PROP_TYPE(TypeTag, PTAG(Problem)) Problem;
+ typedef GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
typedef Dune::FieldVector<Scalar, Grid::dimensionworld> GlobalPosition;
static const int dim = Grid::dimension;
@@ -96,16 +97,15 @@ int main(int argc, char** argv)
////////////////////////////////////////////////////////////
// instantiate and run the concrete problem
////////////////////////////////////////////////////////////
- Dune::FieldVector<double ,dim> L(0);
- Dune::FieldVector<double,dim> H(1);
- Problem problem(gridPtr->leafView(), L, H);
+
+ TimeManager timeManager;
+ Problem problem(timeManager, gridPtr->leafView());
// load restart file if necessarry
if (restart)
problem.deserialize(restartTime);
-
- problem.timeManager().init(problem, 0, dt, tEnd, !restart);
- problem.timeManager().run();
+ timeManager.init(problem, 0, dt, tEnd, !restart);
+ timeManager.run();
return 0;
}
diff --git a/test/decoupled/2p/test_transport_problem.hh b/test/decoupled/2p/test_transport_problem.hh
index ed307bc6bf3c939385f8426aae68a09439566af4..54587edf89610930564d0769c8ccf2a91d0661be 100644
--- a/test/decoupled/2p/test_transport_problem.hh
+++ b/test/decoupled/2p/test_transport_problem.hh
@@ -49,7 +49,7 @@ class TestTransportProblem;
//////////
namespace Properties
{
-NEW_TYPE_TAG(TransportTestProblem, INHERITS_FROM(DecoupledModel, Transport));
+NEW_TYPE_TAG(TransportTestProblem, INHERITS_FROM(DecoupledTwoP, Transport));
// Set the grid type
@@ -123,10 +123,9 @@ SET_SCALAR_PROP(TransportTestProblem, CFLFactor, 1.0);
* where the argument defines the simulation endtime.
*/
template<class TypeTag = TTAG(TransportTestProblem)>
-class TestTransportProblem: public TransportProblem2P<TypeTag, TestTransportProblem<TypeTag> >
+class TestTransportProblem: public TransportProblem2P<TypeTag>
{
- typedef TestTransportProblem<TypeTag> ThisType;
- typedef TransportProblem2P<TypeTag, ThisType> ParentType;
+ typedef TransportProblem2P<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
@@ -134,6 +133,10 @@ class TestTransportProblem: public TransportProblem2P<TypeTag, TestTransportProb
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidState)) FluidState;
+ typedef typename GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
+ typedef typename GET_PROP_TYPE(TypeTag, PTAG(BoundaryTypes)) BoundaryTypes;
+ typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes))::PrimaryVariables PrimaryVariables;
+
enum
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
@@ -141,7 +144,7 @@ class TestTransportProblem: public TransportProblem2P<TypeTag, TestTransportProb
enum
{
- wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
+ wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx, eqIdxSat = Indices::saturationEq
};
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
@@ -149,11 +152,10 @@ class TestTransportProblem: public TransportProblem2P<TypeTag, TestTransportProb
typedef typename GridView::Traits::template Codim<0>::Entity Element;
typedef typename GridView::Intersection Intersection;
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
- typedef Dune::FieldVector<Scalar, dim> LocalPosition;
public:
- TestTransportProblem(const GridView &gridView, const GlobalPosition lowerLeft = 0, const GlobalPosition upperRight = 0) :
- ParentType(gridView), lowerLeft_(lowerLeft), upperRight_(upperRight)
+ TestTransportProblem(TimeManager &timeManager, const GridView &gridView) :
+ ParentType(timeManager, gridView)
{
GlobalPosition vel(0);
vel[0] = 1e-5;
@@ -186,7 +188,7 @@ public:
*
* This problem assumes a temperature of 10 degrees Celsius.
*/
- Scalar temperature(const GlobalPosition& globalPos, const Element& element) const
+ Scalar temperatureAtPos(const GlobalPosition& globalPos) const
{
return 273.15 + 10; // -> 10°C
}
@@ -194,48 +196,64 @@ public:
// \}
- Scalar referencePressure(const GlobalPosition& globalPos, const Element& element) const
+ //! Returns the reference pressure for evaluation of constitutive relations
+ Scalar referencePressureAtPos(const GlobalPosition& globalPos) const
{
return 1e5; // -> 10°C
}
- std::vector<Scalar> source(const GlobalPosition& globalPos, const Element& element)
+ void sourceAtPos(PrimaryVariables &values,const GlobalPosition& globalPos) const
{
- return std::vector<Scalar>(2, 0.0);
+ values = 0;
}
- BoundaryConditions::Flags bctypeSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ /*!
+ * \brief Returns the type of boundary condition.
+ *
+ *
+ * BC for saturation equation can be dirichlet (saturation), neumann (flux), or outflow.
+ */
+ void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos) const
{
- if (globalPos[0] < eps_)
- return Dumux::BoundaryConditions::dirichlet;
- else if (globalPos[0] > upperRight_[0] - eps_)
- return Dumux::BoundaryConditions::outflow;
- else
- return Dumux::BoundaryConditions::neumann;
+ if (globalPos[0] < eps_)
+ {
+ bcTypes.setAllDirichlet();
+ }
+ else if (globalPos[0] > this->bboxMax()[0] - eps_)
+ {
+ bcTypes.setAllOutflow();
+ }
+ // all other boundaries
+ else
+ {
+ bcTypes.setAllNeumann();
+ }
}
- Scalar dirichletSat(const GlobalPosition& globalPos, const Intersection& intersection) const
+ //! set dirichlet condition (saturation [-])
+ void dirichletAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
+ values = 0;
if (globalPos[0] < eps_)
- return 1.0;
- // all other boundaries
- return 0.0;
+ {
+ values = 1.0;
+ }
}
- std::vector<Scalar> neumann(const GlobalPosition& globalPos, const Intersection& intersection) const
+ //! set neumann condition for phases (flux, [kg/(m^2 s)])
+ void neumannAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
- return std::vector<Scalar>(2, 0.0);
+ values = 0;
}
- Scalar initSat(const GlobalPosition& globalPos, const Element& element) const
+ //! return initial solution
+ void initialAtPos(PrimaryVariables &values,
+ const GlobalPosition &globalPos) const
{
- return 0.0;
+ values = 0;
}
private:
- GlobalPosition lowerLeft_;
- GlobalPosition upperRight_;
-
static constexpr Scalar eps_ = 1e-6;
};
} //end namespace