Commit c2298f7f authored by Timo Koch's avatar Timo Koch
Browse files

Merge branch 'feature/improve-fluidsystem-ex' into 'master'

Feature/improve fluidsystem ex

Closes dumux#513

See merge request !4
parents 8efdd356 97d8ac45
......@@ -19,15 +19,15 @@
/*!
* \file
*
* \brief Tutorial problem for a fully coupled twophase box model.
* \brief Tutorial problem for a fully coupled two phase-two component box model.
*/
#ifndef DUMUX_EXERCISE_THREE_B_PROBLEM_HH
#define DUMUX_EXERCISE_THREE_B_PROBLEM_HH
#ifndef DUMUX_EXERCISE_FLUIDSYSTEM_B_PROBLEM_HH
#define DUMUX_EXERCISE_FLUIDSYSTEM_B_PROBLEM_HH
// The numerical model
#include <dumux/porousmediumflow/2p2c/model.hh>
//The box discretization
// The box discretization
#include <dumux/discretization/box/properties.hh>
// The base porous media box problem
......@@ -41,32 +41,31 @@
namespace Dumux{
// Forward declaration of the problem class
template <class TypeTag> class ExerciseThreeProblemTwoPTwoC;
template <class TypeTag> class ExerciseFluidsystemProblemTwoPTwoC;
namespace Properties {
// Create a new type tag for the problem
NEW_TYPE_TAG(ExerciseThreeTwoPTwoCTypeTag, INHERITS_FROM(TwoPTwoC));
NEW_TYPE_TAG(ExerciseThreeBoxTwoPTwoCTypeTag, INHERITS_FROM(BoxModel, ExerciseThreeTwoPTwoCTypeTag));
NEW_TYPE_TAG(ExerciseFluidsystemTwoPTwoCTypeTag, INHERITS_FROM(TwoPTwoC, BoxModel));
// Set the "Problem" property
SET_TYPE_PROP(ExerciseThreeTwoPTwoCTypeTag, Problem, ExerciseThreeProblemTwoPTwoC<TypeTag>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTwoCTypeTag, Problem, ExerciseFluidsystemProblemTwoPTwoC<TypeTag>);
// Set the spatial parameters
SET_TYPE_PROP(ExerciseThreeTwoPTwoCTypeTag, SpatialParams,
ExerciseThreeSpatialParams<typename GET_PROP_TYPE(TypeTag, FVGridGeometry),
SET_TYPE_PROP(ExerciseFluidsystemTwoPTwoCTypeTag, SpatialParams,
ExerciseFluidsystemSpatialParams<typename GET_PROP_TYPE(TypeTag, FVGridGeometry),
typename GET_PROP_TYPE(TypeTag, Scalar)>);
// Set grid and the grid creator to be used
#if HAVE_DUNE_ALUGRID
SET_TYPE_PROP(ExerciseThreeTwoPTwoCTypeTag, Grid, Dune::ALUGrid</*dim=*/2, 2, Dune::cube, Dune::nonconforming>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTwoCTypeTag, Grid, Dune::ALUGrid</*dim=*/2, 2, Dune::cube, Dune::nonconforming>);
#elif HAVE_UG
SET_TYPE_PROP(ExerciseThreeTwoPTwoCTypeTag, Grid, Dune::UGGrid<2>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTwoCTypeTag, Grid, Dune::UGGrid<2>);
#else
SET_TYPE_PROP(ExerciseThreeTwoPTwoCTypeTag, Grid, Dune::YaspGrid<2>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTwoCTypeTag, Grid, Dune::YaspGrid<2>);
#endif // HAVE_DUNE_ALUGRID
// The fluid system property
SET_PROP(ExerciseThreeTwoPTwoCTypeTag, FluidSystem)
SET_PROP(ExerciseFluidsystemTwoPTwoCTypeTag, FluidSystem)
{
private:
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
......@@ -79,10 +78,10 @@ public:
/*!
* \ingroup TwoPBoxModel
*
* \brief Tutorial problem for a fully coupled twophase box model.
* \brief Tutorial problem for a fully coupled two phase-two component box model.
*/
template <class TypeTag>
class ExerciseThreeProblemTwoPTwoC : public PorousMediumFlowProblem<TypeTag>
class ExerciseFluidsystemProblemTwoPTwoC : public PorousMediumFlowProblem<TypeTag>
{
using ParentType = PorousMediumFlowProblem<TypeTag>;
using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
......@@ -104,12 +103,12 @@ class ExerciseThreeProblemTwoPTwoC : public PorousMediumFlowProblem<TypeTag>
using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
public:
ExerciseThreeProblemTwoPTwoC(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
ExerciseFluidsystemProblemTwoPTwoC(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
: ParentType(fvGridGeometry)
, eps_(3e-6)
{
#if !(HAVE_DUNE_ALUGRID || HAVE_UG)
std::cout << "If you want to use simplices instead of cubes, install and use dune-ALUGrid or UGGrid." << std::endl;
std::cout << "If you want to use simplices instead of cubes, install and use dune-ALUGrid or UGGrid." << std::endl;
#endif // !(HAVE_DUNE_ALUGRID || HAVE_UG)
// initialize the fluid system
......@@ -117,9 +116,6 @@ public:
// set the depth of the bottom of the reservoir
depthBOR_ = this->fvGridGeometry().bBoxMax()[dimWorld-1];
// name of the problem and output file
name_ = getParam<std::string>("Problem.Name");
}
/*!
......@@ -127,15 +123,6 @@ public:
*/
// \{
/*!
* \brief Returns the problem name
*
* This is used as a prefix for files generated by the simulation.
*/
const std::string name() const
{ return name_; }
/*!
* \brief Returns the temperature \f$ K \f$
*/
......@@ -174,11 +161,8 @@ public:
*/
PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
{
PrimaryVariables priVars;
priVars.setState(Indices::firstPhaseOnly);
priVars[Indices::pressureIdx] = 200.0e3 + 9.81*1000*(depthBOR_ - globalPos[dimWorld-1]);
priVars[Indices::switchIdx] = 0.0; // 0 % oil saturation on left boundary
return priVars;
// use initial as Dirichlet conditions
return initialAtPos(globalPos);
}
/*!
......@@ -196,14 +180,17 @@ public:
{
// initialize values to zero, i.e. no-flow Neumann boundary conditions
PrimaryVariables values(0.0);
Scalar up = this->fvGridGeometry().bBoxMax()[dimWorld-1];
// extraction of oil on the right boundary for approx. 1.e6 seconds
if (globalPos[dimWorld-1] > up - eps_ && globalPos[0] > 20 && globalPos[0] < 40) {
// oil outflux of 30 g/(m * s) on the right boundary.
// extraction of oil (30 g/m/s) on a segment of the upper boundary
if (globalPos[dimWorld-1] > up - eps_ && globalPos[0] > 20 && globalPos[0] < 40)
{
// we solve for the mole balance, so we have to divide by the molar mass
values[Indices::conti0EqIdx + FluidSystem::H2OIdx] = 0;
values[Indices::conti0EqIdx + FluidSystem::NAPLIdx] = -3e-2/FluidSystem::MyCompressibleComponent::molarMass();
} else {
}
else
{
// no-flow on the remaining Neumann-boundaries.
values[Indices::conti0EqIdx + FluidSystem::H2OIdx] = 0;
values[Indices::conti0EqIdx + FluidSystem::NAPLIdx] = 0;
......@@ -231,8 +218,13 @@ public:
{
PrimaryVariables values(0.0);
// tell the primary variables the phase state, i.e. which phase/phases
// is/are present, because this changes the meaning of the primary variable
// value at the index Indices::switchIdx
values.setState(Indices::firstPhaseOnly);
// use hydrostatic pressure distribution with 2 bar at the top and zero saturation
values[Indices::pressureIdx] = 200.0e3 + 9.81*1000*(depthBOR_ - globalPos[dimWorld-1]); // 200 kPa = 2 bar
values[Indices::switchIdx] = 0.0;
......@@ -250,20 +242,17 @@ public:
*/
PrimaryVariables sourceAtPos(const GlobalPosition &globalPos) const
{
// we do not define any sources
PrimaryVariables values(0.0);
return values;
}
private:
// small epsilon value
Scalar eps_;
// depth at the bottom of the reservoir
Scalar depthBOR_;
std::string name_;
Scalar eps_; //! small epsilon value
Scalar depthBOR_; //! depth at the bottom of the reservoir
};
}
} // end namespace Dumux
#endif
......@@ -21,16 +21,16 @@
*
* \brief Tutorial problem for a fully coupled twophase box model.
*/
#ifndef DUMUX_EXERCISE_THREE_A_PROBLEM_HH
#define DUMUX_EXERCISE_THREE_A_PROBLEM_HH
#ifndef DUMUX_EXERCISE_FLUIDSYSTEM_A_PROBLEM_HH
#define DUMUX_EXERCISE_FLUIDSYSTEM_A_PROBLEM_HH
// The numerical model
#include <dumux/porousmediumflow/2p/model.hh>
//The box discretization
// The box discretization
#include <dumux/discretization/box/properties.hh>
//The grid managers
// The grid managers
#if HAVE_DUNE_ALUGRID
#include <dune/alugrid/grid.hh>
#elif HAVE_UG
......@@ -39,7 +39,7 @@
#include <dune/grid/yaspgrid.hh>
#endif
// The base porous media box problem
// The porous media base problem
#include <dumux/porousmediumflow/problem.hh>
// Spatially dependent parameters
......@@ -61,57 +61,55 @@
namespace Dumux{
// Forward declaration of the problem class
template <class TypeTag> class ExerciseThreeProblemTwoP;
template <class TypeTag> class ExerciseFluidsystemProblemTwoP;
namespace Properties {
// Create a new type tag for the problem
NEW_TYPE_TAG(ExerciseThreeTwoPTypeTag, INHERITS_FROM(TwoP));
NEW_TYPE_TAG(ExerciseThreeBoxTwoPTypeTag, INHERITS_FROM(BoxModel, ExerciseThreeTwoPTypeTag));
NEW_TYPE_TAG(ExerciseFluidsystemTwoPTypeTag, INHERITS_FROM(TwoP, BoxModel));
// Set the "Problem" property
SET_TYPE_PROP(ExerciseThreeTwoPTypeTag, Problem, ExerciseThreeProblemTwoP<TypeTag>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTypeTag, Problem, ExerciseFluidsystemProblemTwoP<TypeTag>);
// Set the spatial parameters
SET_TYPE_PROP(ExerciseThreeTwoPTypeTag, SpatialParams,
ExerciseThreeSpatialParams<typename GET_PROP_TYPE(TypeTag, FVGridGeometry),
SET_TYPE_PROP(ExerciseFluidsystemTwoPTypeTag, SpatialParams,
ExerciseFluidsystemSpatialParams<typename GET_PROP_TYPE(TypeTag, FVGridGeometry),
typename GET_PROP_TYPE(TypeTag, Scalar)>);
// Set grid and the grid creator to be used
// Set grid to be used
#if HAVE_DUNE_ALUGRID
SET_TYPE_PROP(ExerciseThreeTwoPTypeTag, Grid, Dune::ALUGrid</*dim=*/2, 2, Dune::cube, Dune::nonconforming>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTypeTag, Grid, Dune::ALUGrid</*dim=*/2, 2, Dune::cube, Dune::nonconforming>);
#elif HAVE_UG
SET_TYPE_PROP(ExerciseThreeTwoPTypeTag, Grid, Dune::UGGrid<2>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTypeTag, Grid, Dune::UGGrid<2>);
#else
SET_TYPE_PROP(ExerciseThreeTwoPTypeTag, Grid, Dune::YaspGrid<2>);
SET_TYPE_PROP(ExerciseFluidsystemTwoPTypeTag, Grid, Dune::YaspGrid<2>);
#endif // HAVE_DUNE_ALUGRID
// we use the immiscible fluid system here
SET_PROP(ExerciseThreeTwoPTypeTag, FluidSystem)
SET_PROP(ExerciseFluidsystemTwoPTypeTag, FluidSystem)
{
private:
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using TabulatedH2O = Components::TabulatedComponent<Components::H2O<Scalar>>;
using WettingPhase = typename FluidSystems::OnePLiquid<Scalar, TabulatedH2O>;
using LiquidWaterPhase = typename FluidSystems::OnePLiquid<Scalar, TabulatedH2O>;
/*!
* Uncomment first line and comment second line for using the incompressible component
* Uncomment second line and comment first line for using the compressible component
*/
using NonWettingPhase = typename FluidSystems::OnePLiquid<Scalar, MyIncompressibleComponent<Scalar> >;
// using NonWettingPhase = typename FluidSystems::OnePLiquid<Scalar, MyCompressibleComponent<Scalar> >;
using LiquidMyComponentPhase = typename FluidSystems::OnePLiquid<Scalar, MyIncompressibleComponent<Scalar> >;
// using LiquidMyComponent = typename FluidSystems::OnePLiquid<Scalar, MyCompressibleComponent<Scalar> >;
public:
using type = typename FluidSystems::TwoPImmiscible<Scalar, WettingPhase, NonWettingPhase>;
using type = typename FluidSystems::TwoPImmiscible<Scalar, LiquidWaterPhase, LiquidMyComponentPhase>;
};
}
/*!
* \ingroup TwoPBoxModel
*
* \brief Tutorial problem for a fully coupled twophase box model.
*/
template <class TypeTag>
class ExerciseThreeProblemTwoP : public PorousMediumFlowProblem<TypeTag>
class ExerciseFluidsystemProblemTwoP : public PorousMediumFlowProblem<TypeTag>
{
using ParentType = PorousMediumFlowProblem<TypeTag>;
using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
......@@ -140,28 +138,24 @@ class ExerciseThreeProblemTwoP : public PorousMediumFlowProblem<TypeTag>
};
public:
ExerciseThreeProblemTwoP(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
ExerciseFluidsystemProblemTwoP(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
: ParentType(fvGridGeometry)
, eps_(3e-6)
, eps_(3e-6)
{
#if !(HAVE_DUNE_ALUGRID || HAVE_UG)
std::cout << "If you want to use simplices instead of cubes, install and use dune-ALUGrid or UGGrid." << std::endl;
std::cout << "If you want to use simplices instead of cubes, install and use dune-ALUGrid or UGGrid." << std::endl;
#endif // !(HAVE_DUNE_ALUGRID || HAVE_UG)
// initialize the tables for the water properties
std::cout << "Initializing the tables for the water properties" << std::endl;
Components::TabulatedComponent<Components::H2O<Scalar>>::init(/*tempMin=*/273.15,
/*tempMax=*/623.15,
/*numTemp=*/100,
/*pMin=*/0.0,
/*pMax=*/20e6,
/*numP=*/200);
// set the depth of the bottom of the reservoir
depthBOR_ = this->fvGridGeometry().bBoxMax()[dimWorld-1];
// initialize the tables for the water properties
Components::TabulatedComponent<Components::H2O<Scalar>>::init(/*tempMin=*/273.15,
/*tempMax=*/623.15,
/*numTemp=*/100,
/*pMin=*/0.0,
/*pMax=*/20e6,
/*numP=*/200);
// name of the problem and output file
name_ = getParam<std::string>("Problem.Name");
// set the depth of the bottom of the reservoir
depthBOR_ = this->fvGridGeometry().bBoxMax()[dimWorld-1];
}
/*!
......@@ -169,15 +163,6 @@ public:
*/
// \{
/*!
* \brief Returns the problem name
*
* This is used as a prefix for files generated by the simulation.
*/
const std::string name() const
{ return name_; }
/*!
* \brief Returns the temperature \f$ K \f$
*/
......@@ -214,13 +199,10 @@ public:
*
* \param globalPos The global position
*/
PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
PrimaryVariables dirichletAtPos(const GlobalPosition& globalPos) const
{
PrimaryVariables priVars;
priVars[waterPressureIdx] = 200.0e3 + 9.81*1000*(depthBOR_ - globalPos[dimWorld-1]); // 200 kPa = 2 bar
priVars[naplSaturationIdx] = 0.0; // 0 % oil saturation on left boundary
return priVars;
// use the initial values as Dirichlet values
return initialAtPos(globalPos);
}
/*!
......@@ -234,18 +216,21 @@ public:
* For this method, the \a values parameter stores the mass flux
* in normal direction of each phase. Negative values mean influx.
*/
PrimaryVariables neumannAtPos(const GlobalPosition &globalPos) const
PrimaryVariables neumannAtPos(const GlobalPosition& globalPos) const
{
// initialize values to zero, i.e. no-flow Neumann boundary conditions
PrimaryVariables values(0.0);
Scalar up = this->fvGridGeometry().bBoxMax()[dimWorld-1];
// extraction of oil on the right boundary for approx. 1.e6 seconds
if (globalPos[dimWorld-1] > up - eps_ && globalPos[0] > 20 && globalPos[0] < 40) {
// oil outflux of 30 g/(m * s) on the right boundary.
// influx of oil (30 g/m/s) over a segment of the top boundary
if (globalPos[dimWorld-1] > up - eps_ && globalPos[0] > 20 && globalPos[0] < 40)
{
values[contiWEqIdx] = 0;
values[contiNEqIdx] = -3e-2;
} else {
}
else
{
// no-flow on the remaining Neumann-boundaries.
values[contiWEqIdx] = 0;
values[contiNEqIdx] = 0;
......@@ -269,13 +254,14 @@ public:
* For this method, the \a values parameter stores primary
* variables.
*/
PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const
{
PrimaryVariables values(0.0);
values[waterPressureIdx] = 200.0e3 + 9.81*1000*(depthBOR_ - globalPos[dimWorld-1]); // 200 kPa = 2 bar (pw)
values[naplSaturationIdx] = 0.0; // (sn)
// use hydrostatic pressure distribution with 2 bar at the top and zero saturation
values[waterPressureIdx] = 200.0e3 + 9.81*1000*(depthBOR_ - globalPos[dimWorld-1]);
values[naplSaturationIdx] = 0.0;
return values;
}
......@@ -288,20 +274,16 @@ public:
* \f$ [ \textnormal{unit of primary variable} / (m^\textrm{dim} \cdot s )] \f$
* \param globalPos The global position
*/
PrimaryVariables sourceAtPos(const GlobalPosition &globalPos) const
PrimaryVariables sourceAtPos(const GlobalPosition& globalPos) const
{
// we define no source term here
PrimaryVariables values(0.0);
return values;
}
private:
// small epsilon value
Scalar eps_;
std::string name_; //! Problem name
// depth at the bottom of the reservoir
Scalar depthBOR_;
Scalar eps_; //! small epsilon value
Scalar depthBOR_; //! depth at the bottom of the reservoir
};
}
......
# executables for exercise part a & b
#part a: 2pproblem
dune_add_test(NAME exercise3_a
SOURCES exercise3.cc
COMPILE_DEFINITIONS TYPETAG=ExerciseThreeBoxTwoPTypeTag
dune_add_test(NAME exercise-fluidsystem_a
SOURCES exercise-fluidsystem.cc
COMPILE_DEFINITIONS TYPETAG=ExerciseFluidsystemTwoPTypeTag
COMPILE_ONLY # for testing purposes, ignore for the exercise
CMD_ARGS exercise3_a.input)
CMD_ARGS exercise-fluidsystem_a.input)
#part b: 2p2cproblem
dune_add_test(NAME exercise3_b
SOURCES exercise3.cc
COMPILE_DEFINITIONS TYPETAG=ExerciseThreeBoxTwoPTwoCTypeTag
dune_add_test(NAME exercise-fluidsystem_b
SOURCES exercise-fluidsystem.cc
COMPILE_DEFINITIONS TYPETAG=ExerciseFluidsystemTwoPTwoCTypeTag
COMPILE_ONLY # for testing purposes, ignore for the exercise
CMD_ARGS exercise3_b.input)
CMD_ARGS exercise-fluidsystem_b.input)
# add tutorial to the common target
add_dependencies(test_exercises exercise3_a exercise3_b)
# add exercises to the common target
add_dependencies(test_exercises exercise-fluidsystem_a exercise-fluidsystem_b)
# add symlinks for the input files
add_input_file_links()
# Exercise #3 (DuMuX course)
# Exercise Fluidsystem
The aim of this exercise is to get familiar with the _DuMuX_ way of implementing new components (fluids) and fluid systems (mixtures). In the scope of this exercise, a new fictitious component is implemented (exercise _3a_) as well as its mixture with water (exercise _3b_).
The aim of this exercise is to get familiar with the _DuMuX_ way of implementing new components (fluids) and fluid systems (mixtures). In the scope of this exercise, a new fictitious component is implemented
(exercise-fluidsystem a) as well as its mixture with water (exercise-fluidsystem b).
## Problem set-up
The domain has a size of 60 x 60 m and contains two low-permeable lenses. Initially, the domain is fully water saturated and the fictitious component is injected through the middle portion of the upper boundary by means of a Neumann boundary condition. The remaining parts of the upper and the entire lower boundary are Neumann no-flow while on the two lateral sides Dirichlet boundary conditions are applied (hydrostatic conditions for the pressure and zero saturation).
![](../extradoc/exercise3_setup.png)
![](../extradoc/exercise-fluidsystem_setup.png)
## Preparing the exercise
* Navigate to the directory `dumux/tutorial/ex3`
* Navigate to the directory `dumux-course/exercises/exercise-fluidsystem`
### 1. Getting familiar with the code
Locate all the files you will need for this exercise
* The shared __main file__ : `exercise3.cc`
* The __input file__ for part a: `exercise3_a.input`
* The shared __main file__ : `exercise-fluidsystem.cc`
* The __input file__ for part a: `exercise-fluidsystem_a.input`
* The __problem file__ for part a: `2pproblem.hh`
* The __input file__ for part b: `exercise3_b.input`
* The __input file__ for part b: `exercise-fluidsystem_b.input`
* The __problem file__ for part b: `2p2cproblem.hh`
* The __spatial parameters file__: `spatialparams.hh`
......@@ -34,34 +35,39 @@ To see more components, fluidsystems and binarycoefficients implementations, hav
In the following, the basic steps required to set the desired fluid system are outlined. Here, this is done in the __problem file__, i.e. for this part of the exercise the code shown below is taken from the `2pproblem.hh` file.
In this part of the exercise we will consider a system consisting of two immiscible phases. Therefore, the _TypeTag_ for this problem (`ExerciseThreeBoxTwoPTypeTag`) derives from a base _TypeTag_ (`ExerciseThreeTwoPTypeTag`) that itself derives from the `TwoP` _TypeTag_ (immiscible two-phase model properties).
In this part of the exercise we will consider a system consisting of two immiscible phases. Therefore, the _TypeTag_ for this problem (`ExerciseFluidsystemTwoPTypeTag`) derives from
the `TwoP` _TypeTag_ (immiscible two-phase model properties) and the `BoxModel` _TypeTag_ (specifies properties of the discretization scheme).
```c++
NEW_TYPE_TAG(ExerciseThreeTwoPTypeTag, INHERITS_FROM(TwoP));
NEW_TYPE_TAG(ExerciseFluidsystemTwoPTypeTag, INHERITS_FROM(TwoP, BoxModel));
```
In order to be able to derive from this _TypeTag_, the declaration of the `TwoP` _TypeTag_ has to be included. It can be found in the `2p/model.hh` header:
In order to be able to derive from these _TypeTags_, the declarations of the `TwoP` _TypeTag_ and `BoxModel` _TypeTag_ have to be included.
The `TwoP` _TypeTag_ can be found in the `2p/model.hh` header:
```c++
// The numerical model
#include <dumux/porousmediumflow/2p/model.hh>
```
Additionally, the _TypeTag_ for this problem (`ExerciseThreeBoxTwoPTypeTag`) derives from the `BoxModel` _TypeTag_, to specify properties of the discretization scheme. For a cell-centered scheme, you could derive from `CCTpfaModel` or `CCMpfaModel` instead. Again the corresponding header has to be included
while the `BoxModel` _TypeTag_ can be found in the `box/properties.hh` header:
```c++
// The discretization
#include <dumux/discretization/box/properties.hh>
```
As wetting phase we want to use water and we want to precompute tables on which the properties are then interpolated in order to save computational time. Thus, in a first step we have to include the following headers:
For a cell-centered scheme, you could derive from `CCTpfaModel` or `CCMpfaModel` instead (and, of course, include the right headers).
As one of the two phases we want to use water and we want to precompute tables on which the properties are then interpolated in order to save computational time. Thus, in a first step we have to include the following headers:
```c++
// The water component
#include <dumux/material/components/tabulatedcomponent.hh>
#include <dumux/material/components/h2o.hh>
```
The non-wetting phase will be our new component, where we want to implement an incompressible and a compressible variant. The respective headers are prepared, but still incomplete. The compressible variant is still commented so that compilation does not fail when finishing the incompressible variant.
The other phase will be created from our new component, where we want to implement an incompressible and a compressible variant.
The respective headers are prepared, but still incomplete. The compressible variant is still commented so that compilation does not fail when finishing the incompressible variant.
```c++
// The components that will be created in this exercise
......@@ -82,26 +88,28 @@ This fluid system expects __phases__ as input and so far we have only included t
#include <dumux/material/fluidsystems/1pliquid.hh>
```
which creates a _liquid phase_ from a given component. Finally, using all of the included classes we set the fluid system property by choosing that the non-wetting phase is a one-phase liquid (OnePLiquid) consisting of the incompressible fictitious component and that the wetting-phase consists of tabulated water in the immiscible fluid system:
which creates a _liquid phase_ from a given component. Finally, using all of the included classes we set the fluid system property by choosing that
the water phase is liquid (OnePLiquid) and consists of the tabulated water component, and
the other phase is liquid as well and consists of the incompressible fictitious component. Both will make up the immiscible fluid system (TwoPImmiscible):
```c++
// we use the immiscible fluid system here
SET_PROP(ExerciseThreeTwoPTypeTag, FluidSystem)
SET_PROP(ExerciseFluidsystemTwoPTypeTag, FluidSystem)
{
private:
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using TabulatedH2O = Components::TabulatedComponent<Components::H2O<Scalar>>;
using WettingPhase = typename FluidSystems::OnePLiquid<Scalar, TabulatedH2O>;
using LiquidWaterPhase = typename FluidSystems::OnePLiquid<Scalar, TabulatedH2O>;
/*!
* Uncomment first line and comment second line for using the incompressible component
* Uncomment second line and comment first line for using the compressible component
*/
using NonWettingPhase = typename FluidSystems::OnePLiquid<Scalar, MyIncompressibleComponent<Scalar> >;
// using NonWettingPhase = typename FluidSystems::OnePLiquid<Scalar, MyCompressibleComponent<Scalar> >;
using LiquidMyComponentPhase = typename FluidSystems::OnePLiquid<Scalar, MyIncompressibleComponent<Scalar> >;
// using LiquidMyComponent = typename FluidSystems::OnePLiquid<Scalar, MyCompressibleComponent<Scalar> >;
public:
using type = typename FluidSystems::TwoPImmiscible<Scalar, WettingPhase, NonWettingPhase>;
using type = typename FluidSystems::TwoPImmiscible<Scalar, LiquidWaterPhase, LiquidMyComponentPhase>;
};
```
......@@ -128,8 +136,8 @@ Implement an incompressible component into the file `myincompressiblecomponent.h
| Parameter | unit | value |
| -----| --------| -------- |
| $`M`$ | $`Kg/mol`$ | $`131.39 \cdot 10^{-3}`$ |
| $`\rho_{liquid}`$ | $`Kg/m^3`$ | $`1460`$ |
| $`M`$ | $`kg/mol`$ | $`131.39 \cdot 10^{-3}`$ |
| $`\rho_{liquid}`$ | $`kg/m^3`$ | $`1460`$ |
| $`\mu_{liquid}`$ | $`Pa \cdot s`$ | $`5.7 \cdot 10^{-4}`$ |
In order to do so, have a look at the files `dumux/material/components/base.hh` and `dumux/material/components/liquid.hh` to see how the interfaces are defined and overload them accordingly.
......@@ -137,14 +145,14 @@ In order to do so, have a look at the files `dumux/material/components/base.hh`
In order to execute the program, change to the build directory and compile and execute the program by typing
```bash
cd build-cmake/tutorial/ex3
make exercise3_a
./exercise3_a exercise3_a.input
cd build-cmake/exercises/exercise-fluidsystem
make exercise-fluidsystem_a
./exercise-fluidsystem_a exercise-fluidsystem_a.input
```
The saturation distribution of the nonwetting phase at the final simulation time should look like this:
The saturation distribution of the nonwetting phase S$`_n`$ (the phase consisting of our fictitious incompressible component) at the final simulation time should look like this:
![](../extradoc/exercise3_a_solution.png)
![](../extradoc/exercise-fluidsystem_a_solution.png)