Commit aae93b02 authored by Benjamin Faigle's avatar Benjamin Faigle
Browse files

Correct and remove decoupled 2p2c doxygen errors.

Improved docu on adaptive models: Everything should be properly doxygened now.
Added Adaptive2p2c doxygen module. 

git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@9203 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent 04dd27c7
......@@ -261,12 +261,16 @@
*/
/*!
* \ingroup IMPET
* \defgroup IMPEC Miscible IMPEC
* \defgroup IMPEC Miscible (Compositional) IMPEC
*/
/*!
* \ingroup IMPEC
* \defgroup multiphase Multiphase Compositional Models
*/
/*!
* \ingroup multiphase
* \defgroup Adaptive2p2c (Grid-)Adaptive Multiphase Compositional Models
*/
/*!
* \ingroup IMPEC
* \defgroup multiphysics Multiphysics Compositional Models
......
......@@ -19,6 +19,7 @@
/*!
* \ingroup IMPEC
* \ingroup Adaptive2p2c
* \ingroup IMPETProperties
*
* \file
......
......@@ -30,7 +30,7 @@
namespace Dumux
{
/*!
* \ingroup IMPET
* \ingroup Adaptive2p2c
*/
//! Class including the data of a grid cell needed if an adaptive grid is used.
/*! The class provides model-specific functions needed to adapt the stored cell data to a new (adapted) grid.
......@@ -73,7 +73,7 @@ private:
public:
//! A container for all necessary variables to map an old solution to a new grid
/*
/**
* If the primary variables (pressure, total concentrations) are mapped to a new grid,
* the secondary variables can be calulated. For the mapping between sons and father, it
* is in addition necessary to know about how many sons live in each father ("count").
......@@ -85,14 +85,15 @@ public:
*/
struct AdaptedValues
{
Dune::FieldVector<Scalar,2> totalConcentration_;
Dune::FieldVector<Scalar,2> pressure_;
Dune::FieldVector<Scalar,2> totalConcentration_; //!< Transport primary variables
Dune::FieldVector<Scalar,2> pressure_; //!< Pressure primary variables
//! Storage for volume derivatives, as transport estimate on old pressure field is not correct after refinement
Dune::FieldVector<Scalar,3> volumeDerivatives_;
Scalar cellVolume;
FluxData2P2C<TypeTag> fluxData_;
int subdomain;
int count;
int isRefined;
Scalar cellVolume; //!< Cell volume for transformation of volume-specific primary variables
FluxData2P2C<TypeTag> fluxData_; //!< Information on old flux direction
int subdomain; //!< subdomain
int count; //!< counts the number of cells averaged
int isRefined; //!< Indicates if cell is refined
AdaptedValues()
{
totalConcentration_=0.;
......
......@@ -181,8 +181,8 @@ public:
protected:
Problem& problem_;
Scalar maxError_; //!> Maximum volume error of all cells
bool enableVolumeIntegral; //!> Enables the volume integral of the pressure equation
Scalar maxError_; //!< Maximum volume error of all cells
bool enableVolumeIntegral; //!< Enables the volume integral of the pressure equation
Scalar ErrorTermFactor_; //!< Handling of error term: relaxation factor
Scalar ErrorTermLowerBound_; //!< Handling of error term: lower bound for error dampening
Scalar ErrorTermUpperBound_; //!< Handling of error term: upper bound for error dampening
......
......@@ -55,7 +55,7 @@ NEW_PROP_TAG(GridAdaptEnableSecondHalfEdge);
}
//! The finite volume model for the solution of the compositional pressure equation
/*! \ingroup multiphase
/*! \ingroup Adaptive2p2c
* Provides a Finite Volume implementation for the pressure equation of a compressible
* system with two components. An IMPES-like method is used for the sequential
* solution of the problem. Diffusion is neglected, capillarity can be regarded.
......@@ -224,9 +224,9 @@ protected:
//! Matrix for vector rotation used in mpfa
DimMatrix R_;
bool enableVolumeIntegral; //!> Enables the volume integral of the pressure equation
bool enableMPFA; //!> Enables mpfa method to calculate the fluxes near hanging nodes
bool enableSecondHalfEdge; //!> If possible, 2 interaction volumes are used for the mpfa method near hanging nodes
bool enableVolumeIntegral; //!< Enables the volume integral of the pressure equation
bool enableMPFA; //!< Enables mpfa method to calculate the fluxes near hanging nodes
bool enableSecondHalfEdge; //!< If possible, 2 interaction volumes are used for the mpfa method near hanging nodes
//! The 2p Mpfa pressure module, that is only used for the calulation of transmissibility of the second interaction volumes
FVMPFAL2PFABoundPressure2PAdaptive<TypeTag>* pressureModelAdaptive2p_;
};
......
......@@ -114,7 +114,8 @@ class FVPressure2P2CMultiPhysics : public FVPressure2P2C<TypeTag>
typedef Dune::FieldVector<Scalar, GET_PROP_VALUE(TypeTag, NumPhases)> PhaseVector;
typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
//! @copydoc FVPressure::EntryType
typedef Dune::FieldVector<Scalar, 2> EntryType;
//! Access functions to the current problem object
Problem& problem()
{ return this->problem_; }
......@@ -124,13 +125,13 @@ public:
//function which assembles the system of equations to be solved
void assemble(bool first);
void get1pSource(Dune::FieldVector<Scalar, 2>&, const Element&, const CellData&);
void get1pSource(EntryType&, const Element&, const CellData&);
void get1pStorage(Dune::FieldVector<Scalar, 2>&, const Element&, CellData&);
void get1pStorage(EntryType&, const Element&, CellData&);
void get1pFlux(Dune::FieldVector<Scalar, 2>&, const Intersection&, const CellData&);
void get1pFlux(EntryType&, const Intersection&, const CellData&);
void get1pFluxOnBoundary(Dune::FieldVector<Scalar, 2>&,
void get1pFluxOnBoundary(EntryType&,
const Intersection&, const CellData&);
//initialize mult-physics-specific pressure model stuff
......@@ -189,7 +190,7 @@ protected:
Dune::BlockVector<Dune::FieldVector<int,1> > nextSubdomain; //! vector holding next subdomain
const GlobalPosition& gravity_; //!< vector including the gravity constant
static constexpr int pressureType = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w \f$, \f$ 1 = p_n \f$, \f$ 2 = p_{global} \f$)
Dune::Timer timer_; //!> A timer for the time spent on the multiphysics framework.
Dune::Timer timer_; //!< A timer for the time spent on the multiphysics framework.
//! Indices of matrix and rhs entries
/**
......@@ -306,7 +307,6 @@ void FVPressure2P2CMultiPhysics<TypeTag>::assemble(bool first)
* \param sourceEntry The Matrix and RHS entries
* \param elementI The element I
* \param cellDataI Data of cell I
* \param first Flag if pressure field is unknown
*/
template<class TypeTag>
void FVPressure2P2CMultiPhysics<TypeTag>::get1pSource(Dune::FieldVector<Scalar, 2>& sourceEntry, const Element& elementI, const CellData& cellDataI)
......@@ -823,6 +823,7 @@ void FVPressure2P2CMultiPhysics<TypeTag>::updateMaterialLaws(bool postTimeStep)
* primary variables. Only a simple flash calulation has to be carried out,
* as phase distribution is already known: single-phase.
* \param elementI The element
* \param cellData The cell data of the current element
* \param postTimeStep Flag indicating if we have just completed a time step
*/
template<class TypeTag>
......
......@@ -31,7 +31,7 @@
namespace Dumux
{
//! Miscible Transport step in a Finite Volume discretization
//! Compositional transport step in a Finite Volume discretization
/*! \ingroup multiphase
* The finite volume model for the solution of the transport equation for compositional
* two-phase flow.
......@@ -98,7 +98,9 @@ class FVTransport2P2C
typedef Dune::FieldMatrix<Scalar,dim,dim> DimMatrix;
typedef Dune::FieldVector<Scalar, NumPhases> PhaseVector;
typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
protected:
//! @copydoc FVPressure::EntryType
typedef Dune::FieldVector<Scalar, 2> EntryType;
//! Acess function for the current problem
Problem& problem()
{return problem_;};
......@@ -109,11 +111,11 @@ public:
void updateTransportedQuantity(TransportSolutionType& updateVector);
// Function which calculates the flux update
void getFlux(Dune::FieldVector<Scalar, 2>&, Dune::FieldVector<Scalar, 2>&,
void getFlux(EntryType&, EntryType&,
const Intersection&, CellData&);
// Function which calculates the boundary flux update
void getFluxOnBoundary(Dune::FieldVector<Scalar, 2>&, Dune::FieldVector<Scalar, 2>&,
void getFluxOnBoundary(EntryType&, EntryType&,
const Intersection&, const CellData&);
void evalBoundary(GlobalPosition,const Intersection&,FluidState &, PhaseVector &);
......@@ -128,7 +130,7 @@ public:
totalConcentration_[nCompIdx].resize(problem_.gridView().size(0));
};
//! \brief Write data files
//! \brief Write transport variables into the output files
/* \param writer applied VTK-writer */
template<class MultiWriter>
void addOutputVtkFields(MultiWriter &writer)
......@@ -164,7 +166,7 @@ public:
/*! \name Access functions for protected variables */
//@{
//! Return the vector of the transported quantity
/*! For an immiscible IMPES scheme, this is the saturation. For Miscible simulations, however,
/*! For an immiscible IMPES scheme, this is the saturation. For compositional simulations, however,
* the total concentration of all components is transported.
*/
TransportSolutionType& transportedQuantity() DUNE_DEPRECATED
......@@ -189,7 +191,7 @@ public:
//@}
//! Constructs a FVTransport2P2C object
/*!
* Currently, the miscible transport scheme can not be applied with a global pressure / total velocity
* Currently, the compositional transport scheme can not be applied with a global pressure / total velocity
* formulation.
*
* \param problem a problem class object
......@@ -209,13 +211,13 @@ public:
}
protected:
TransportSolutionType totalConcentration_;
TransportSolutionType totalConcentration_; //!< private vector of transported primary variables
Problem& problem_;
bool impet_;
int averagedFaces_;
bool impet_; //!< indicating if we are in an estimate (false) or real impet (true) step.
int averagedFaces_; //!< number of faces were flux was restricted
static const int pressureType = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w \f$, \f$ 1 = p_n \f$, \f$ 2 = p_{global} \f$)
int restrictFluxInTransport_;
int restrictFluxInTransport_; //!< Restriction of flux on new pressure field if direction reverses from the pressure equation
bool switchNormals;
};
......@@ -225,12 +227,11 @@ protected:
* \f[
C^{\kappa , new} = C^{\kappa , old} + u,
* \f]
* where \f$ u = \sum_{element faces} \boldsymbol{v}_{\alpha} * \varrho_{\alpha} * X^{\kappa}_{\alpha} * \boldsymbol{n} * A_{element face} \f$,
* \f$ \boldsymbol{n} \f$ is the face normal and \f$ A_{element face} \f$ is the face area.
* where \f$ u = \sum_{\gamma} \boldsymbol{v}_{\alpha} * \varrho_{\alpha} * X^{\kappa}_{\alpha} * \boldsymbol{n} * A_{\gamma} \f$,
* \f$ \boldsymbol{n} \f$ is the face normal and \f$ A_{\gamma} \f$ is the face area of face \f$ \gamma \f$.
*
* In addition to the \a update vector, the recommended time step size \a dt is calculated
* employing a CFL condition.
* This method = old concentrationUpdate()
*
* \param t Current simulation time \f$\mathrm{[s]}\f$
* \param[out] dt Time step size \f$\mathrm{[s]}\f$
......@@ -340,14 +341,15 @@ void FVTransport2P2C<TypeTag>::updateTransportedQuantity(TransportSolutionType&
}
//! Get flux at an interface between two cells
/** The flux thorugh \f$ \gamma{ij} \f$ is calculated according to the underlying pressure field,
/** The flux through \f$ \gamma \f$ is calculated according to the underlying pressure field,
* calculated by the pressure model.
* \f[ - A_{\gamma_{ij}} \cdot \mathbf{u} \cdot (\mathbf{n}_{\gamma_{ij}} \cdot \mathbf{u})\cdot \mathbf{K}
\sum_{\alpha} \varrho_{\alpha} \lambda_{\alpha} \sum_{\kappa} X^{\kappa}_{\alpha}
\left( \frac{p_{\alpha,j}^t - p^{t}_{\alpha,i}}{\Delta x} + \varrho_{\alpha} \mathbf{g}\right) \f]
* Here, \f$ \mathbf{u} \f$ is the normalized vector connecting the cell centers, and \f$ \mathbf{n}_{\gamma_{ij}} \f$
* represents the normal of the face \f$ \gamma{ij} \f$. Due to the nature of the Primay Variable, the (volume-)specific
* \f[ - A_{\gamma} \mathbf{n}^T_{\gamma} \mathbf{K} \sum_{\alpha} \varrho_{\alpha} \lambda_{\alpha}
\mathbf{d}_{ij} \left( \frac{p_{\alpha,j}^t - p^{t}_{\alpha,i}}{\Delta x} + \varrho_{\alpha} \mathbf{g}^T \mathbf{d}_{ij} \right)
\sum_{\kappa} X^{\kappa}_{\alpha} \f]
* Here, \f$ \mathbf{d}_{ij} \f$ is the normalized vector connecting the cell centers, and \f$ \mathbf{n}_{\gamma} \f$
* represents the normal of the face \f$ \gamma \f$. Due to the nature of the primay Variable, the (volume-)specific
* total mass concentration, this represents a mass flux per cell volume.
*
* \param fluxEntries The flux entries, mass influx from cell \f$j\f$ to \f$i\f$.
* \param timestepFlux flow velocities for timestep estimation
* \param intersection The intersection
......@@ -587,11 +589,11 @@ void FVTransport2P2C<TypeTag>::getFlux(Dune::FieldVector<Scalar, 2>& fluxEntries
return;
}
//! Get flux on Boundary
/** The flux thorugh \f$ \gamma{ij} \f$ is calculated according to the underlying pressure field,
/** The flux through \f$ \gamma \f$ is calculated according to the underlying pressure field,
* calculated by the pressure model.
* \f[ - A_{\gamma_{ij}} \cdot \mathbf{u} \cdot (\mathbf{n}_{\gamma_{ij}} \cdot \mathbf{u})\cdot \mathbf{K}
* \f[ - A_{\gamma} \mathbf{n}^T_{\gamma} \mathbf{K} \mathbf{d}_{i-Boundary}
\sum_{\alpha} \varrho_{\alpha} \lambda_{\alpha} \sum_{\kappa} X^{\kappa}_{\alpha}
\left( \frac{p_{\alpha,j}^t - p^{t}_{\alpha,i}}{\Delta x} + \varrho_{\alpha} \mathbf{g}\right) \f]
\left( \frac{p_{\alpha,Boundary}^t - p^{t}_{\alpha,i}}{\Delta x} + \varrho_{\alpha}\mathbf{g}^T \mathbf{d}_{i-Boundary}\right) \f]
* Here, \f$ \mathbf{u} \f$ is the normalized vector connecting the cell centers, and \f$ \mathbf{n}_{\gamma_{ij}} \f$
* represents the normal of the face \f$ \gamma{ij} \f$.
* \param fluxEntries The flux entries, mass influx from cell \f$j\f$ to \f$i\f$.
......@@ -805,7 +807,7 @@ void FVTransport2P2C<TypeTag>::getFluxOnBoundary(Dune::FieldVector<Scalar, 2>& f
* possible to define boundaries by means of a saturation. This choice determines
* which version of flash calculation is necessary to get to the composition at
* the boundary.
* \param globalPosFace Face of the current boundary
* \param globalPosFace The global position of the face of the current boundary cell
* \param intersection The current intersection
* \param BCfluidState FluidState object that is used for the boundary
* \param pressBound Boundary values of phase pressures
......
......@@ -39,8 +39,8 @@ namespace Properties
// forward declaration of properties
NEW_PROP_TAG(GridAdaptEnableMultiPointFluxApproximation);
}
//! Miscible Transport step in a Finite Volume discretization
/*! \ingroup multiphase
//! Compositional Transport step in a Finite Volume discretization
/*! \ingroup Adaptive2p2c
* The finite volume model for the solution of the transport equation for compositional
* two-phase flow.
* \f[
......@@ -144,12 +144,12 @@ protected:
* \f[
C^{\kappa , new} = C^{\kappa , old} + u,
* \f]
* where \f$ u = \sum_{element faces} \boldsymbol{v}_{\alpha} * \varrho_{\alpha} * X^{\kappa}_{\alpha} * \boldsymbol{n} * A_{element face} \f$,
* \f$ \boldsymbol{n} \f$ is the face normal and \f$ A_{element face} \f$ is the face area.
* where \f$ u = \sum_{\gamma} \boldsymbol{v}_{\alpha} * \varrho_{\alpha} * X^{\kappa}_{\alpha} * \boldsymbol{n} * A_{\gamma} \f$,
* \f$ \boldsymbol{n} \f$ is the face normal and \f$ A_{\gamma} \f$ is the face area of face \f$ \gamma \f$.
*
* In addition to the \a update vector, the recommended time step size \a dt is calculated
* employing a CFL condition.
* This method = old concentrationUpdate()
* employing a CFL condition. This method uses a standard \a Tpfa method for regular fluxes,
* and a \a mpfa can be used near hanging nodes.
*
* \param t Current simulation time \f$\mathrm{[s]}\f$
* \param[out] dt Time step size \f$\mathrm{[s]}\f$
......@@ -261,6 +261,24 @@ void FVTransport2P2CAdaptive<TypeTag>::update(const Scalar t, Scalar& dt, Transp
return;
}
//! Compute flux over an irregular interface using a \a mpfa method
/** A mpfa l-method is applied to calculate fluxes near hanging nodes, using:
* \f[
- \sum_{\alpha} \varrho_{\alpha} \lambda_{\alpha}
\left( \sum_k \tau_{2k} p^t_{\alpha,k} + \varrho_{\alpha} \sum_k \tau_{2k} \mathbf{g}^T \mathbf{x}_{k} \right)
\sum_{\kappa} X^{\kappa}_{\alpha}
\f]
*
* We provide two options: Calculating the flux expressed by twice the flux
* through the one unique interaction region on the hanging node if one
* halfedge is stored (eg on boundaries). Or using the second interaction
* region covering neighboring cells.
*
* \param fluxEntries The flux entries, mass influx from cell \f$j\f$ to \f$i\f$.
* \param timestepFlux flow velocities for timestep estimation
* \param intersectionIterator Iterator of the intersection between cell I and J
* \param cellDataI The cell data for cell \f$i\f$
*/
template<class TypeTag>
void FVTransport2P2CAdaptive<TypeTag>::getMpfaFlux(Dune::FieldVector<Scalar, 2>& fluxEntries, Dune::FieldVector<Scalar, 2>& timestepFlux,
const IntersectionIterator& intersectionIterator, CellData& cellDataI)
......
......@@ -28,7 +28,7 @@
namespace Dumux
{
//! Miscible Transport step in a Finite Volume discretization
//! Compositional Transport Step in a Finite Volume discretization
/*!
* \ingroup multiphysics
* The finite volume model for the solution of the transport equation for compositional
......@@ -40,10 +40,10 @@ namespace Dumux
* \f$ p_{\alpha} \f$ denotes the phase pressure, \f$ \bf{K} \f$ the absolute permeability, \f$ \lambda_{\alpha} \f$ the phase mobility,
* \f$ \rho_{\alpha} \f$ the phase density and \f$ \bf{g} \f$ the gravity constant and \f$ C^{\kappa} \f$ the total Component concentration.
*
* The model domain is automatically divided
* in a single-phase and a two-phase domain. The full 2p2c model is only evaluated within the
* two-phase subdomain, whereas a single-phase transport model is computed in the rest of the
* domain.
* The model domain is automatically divided into a single-phase and a two-phase domain. As the flux computation is relatively cheap,
* the same method is used for the real transport step independently of the subdomain.
* The pressure equation does not need any derivatives in simple
* subdomains, therefore in the transport estimate step inter-cell fluxes in the simple subdomain are omitted.
*
* \tparam TypeTag The Type Tag
*/
......@@ -99,7 +99,6 @@ public:
/**
* \param problem a problem class object
*/
FVTransport2P2CMultiPhysics(Problem& problem) : FVTransport2P2C<TypeTag>(problem)
{}
......@@ -113,12 +112,11 @@ public:
* \f[
C^{\kappa , new} = C^{\kappa , old} + u,
* \f]
* where \f$ u = \sum_{element faces} \boldsymbol{v}_{\alpha} * \varrho_{\alpha} * X^{\kappa}_{\alpha} * \boldsymbol{n} * A_{element face} \f$,
* \f$ \boldsymbol{n} \f$ is the face normal and \f$ A_{element face} \f$ is the face area.
* where \f$ u = \sum_{\gamma} \boldsymbol{v}_{\alpha} * \varrho_{\alpha} * X^{\kappa}_{\alpha} * \boldsymbol{n} * A_{\gamma} \f$,
* \f$ \boldsymbol{n} \f$ is the face normal and \f$ A_{\gamma} \f$ is the face area of face \f$ \gamma \f$.
*
* In addition to the \a update vector, the recommended time step size \a dt is calculated
* employing a CFL condition.
* This method = old concentrationUpdate()
*
* \param t Current simulation time \f$\mathrm{[s]}\f$
* \param[out] dt Time step size \f$\mathrm{[s]}\f$
......
......@@ -36,7 +36,7 @@ namespace Dumux
* It is used in case of a multiphysics approach. For the non-present phase,
* no information is stored but 0-values are returned to allow for general output
* methods.
* The "flash" calulation routines are in the decoupled flash constrain solver, see
* The "flash" calculation routines are in the decoupled flash constrain solver, see
* Dumux::CompositionalFlash .
* \tparam TypeTag The property Type Tag
*/
......@@ -72,6 +72,7 @@ public:
return 0.;
};
//! \brief Returns the index of the phase that is present in that cell.
int presentPhaseIdx() const
{
return presentPhaseIdx_;
......@@ -82,7 +83,11 @@ public:
*/
Scalar pressure(int phaseIdx) const
{ return pressure_[phaseIdx]; }
/*!
* \brief Returns the density of a phase \f$\mathrm{[kg/m^3]}\f$.
*
* \param phaseIdx the index of the phase
*/
Scalar density(int phaseIdx) const
{
if(phaseIdx == presentPhaseIdx_)
......
......@@ -35,7 +35,7 @@
namespace Dumux
{
/*!
* \ingroup IMPET
* \ingroup Adaptive2p2c
*/
//! Class holding additionally mpfa data of adaptive compositional models.
/*!
......@@ -71,51 +71,57 @@ private:
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
};
enum // for first and second half edge 2D / subvolume face 3D
enum //!< for first and second half edge (2D) or subvolume face (3D)
{
first = 0, second = 1
};
// convenience shortcuts for Vectors/Matrices
typedef Dune::FieldVector<Scalar, GridView::dimensionworld> GlobalPosition;
typedef Dune::FieldVector<Scalar,dim+1> T1ransmissivityMatrix;
typedef Dune::FieldVector<Scalar,dim+1> TransmissivityMatrix;
protected:
// in the 2D case, we need to store 1 additional cell. In 3D, we store
// dim-1=2 cells for one interaction volume, and 4 if two interaction volumes
// are regarded.
/** in the 2D case, we need to store 1 additional cell. In 3D, we store
* dim-1=2 cells for one interaction volume, and 4 if two interaction volumes
* are regarded.
*/
const static int storageRequirement = (dim-1)*(dim-1);
//! Storage object for data related to the MPFA method
/*
/**
* This Struct stores the transmissibility coefficients
* for the two half-eges of an irregular interface (one
* near a hanging node) in an h-adaptive simulation.
*/
struct mpfaData
{
T1ransmissivityMatrix T1_[2];
TransmissivityMatrix T1_[2];
int globalIdx3_[storageRequirement];
GlobalPosition globalPos3_[storageRequirement];
std::vector<IntersectionIterator> secondHalfEdgeIntersection_;
bool hasSecondHalfEdge;
//! Constructor for the local storage object of mpfa data
mpfaData()
{
hasSecondHalfEdge = false;
}
//! stores an intersection
/** This also provides the information that both half-edges are
* regarded and information was stored.
* \param is23 Intersection pointing to 3rd cell of additional interaction region
*/
void setIntersection(IntersectionIterator& is23)
{
secondHalfEdgeIntersection_.push_back(is23);
hasSecondHalfEdge = true;
};
//! Acess method to the stored intersection
const IntersectionIterator& getIntersection()
{
return secondHalfEdgeIntersection_[0];
};
};
std::map<IdType, mpfaData> irregularInterfaceMap_;
const Grid& grid_;
std::map<IdType, mpfaData> irregularInterfaceMap_; //!< Storage container for mpfa information
const Grid& grid_; //!< The grid
public:
//! Constructs a VariableClass object
......@@ -153,8 +159,20 @@ public:
problem.pressureModel().adaptPressure();
}
//! Stores Mpfa Data on an intersection
/** The method stores information to the interaction region (Transmissitivity
* as well as details of the 3rd cell in the region) into a storage container.
* The key to each element is the index of the intersection, seen from the smaller
* cell (only this is unique).
* If we arrive from the "wrong" (i.e. non-unique) direction, we invert fluxes.
*
* \param irregularIs The current irregular intersection
* \param T1 Transmissitivity matrix for flux calculations
* \param globalPos3 The position of the 3rd cell of the interaction region
* \param globalIdx3 The index of the 3rd cell of the interaction region
*/
void storeMpfaData(const typename GridView::Intersection & irregularIs,
const T1ransmissivityMatrix& T1,
const TransmissivityMatrix& T1,
const GlobalPosition& globalPos3,
const int& globalIdx3)
{
......@@ -182,11 +200,24 @@ public:
irregularInterfaceMap_[intersectionID].globalIdx3_[0] = globalIdx3;
return;
}
//! Stores Mpfa Data on an intersection for both half-edges
/** The method stores information of both interaction regions (Transmissitivity
* as well as details of the 3rd cells of both regions) into a storage container.
* The key to each element is the index of the intersection, seen from the smaller
* cell (only this is unique).
* If we arrive from the "wrong" (i.e. non-unique) direction, we invert fluxes.
*
* \param irregularIs The current irregular intersection
* \param secondHalfEdgeIntersectionIt Iterator to the intersection connecting the second interaction region
* \param T1 Transmissitivity matrix for flux calculations: unique interaction region
* \param T1_secondHalfEdge Second transmissitivity matrix for flux calculations for non-unique interaction region
* \param globalPos3 The position of the 3rd cell of the first interaction region
* \param globalIdx3 The index of the 3rd cell of the first interaction region
*/
void storeMpfaData(const typename GridView::Intersection & irregularIs,
IntersectionIterator& secondHalfEdgeIntersectionIt,
const T1ransmissivityMatrix& T1,
const T1ransmissivityMatrix& T11_secondHalfEdge,
const TransmissivityMatrix& T1,
const TransmissivityMatrix& T1_secondHalfEdge,
const GlobalPosition& globalPos3,
const int& globalIdx3)
{
......@@ -207,16 +238,16 @@ public:
irregularInterfaceMap_[intersectionID].T1_[first][1] = - T1[1];
irregularInterfaceMap_[intersectionID].T1_[first][0] = - T1[2];
irregularInterfaceMap_[intersectionID].T1_[second][2] = - T11_secondHalfEdge[0];
irregularInterfaceMap_[intersectionID].T1_[second][1] = - T11_secondHalfEdge[1];
irregularInterfaceMap_[intersectionID].T1_[second][0] = - T11_secondHalfEdge[2];
irregularInterfaceMap_[intersectionID].T1_[second][2] = - T1_secondHalfEdge[0];
irregularInterfaceMap_[intersectionID].T1_[second][1] = - T1_secondHalfEdge[1];
irregularInterfaceMap_[intersectionID].T1_[second][0] = - T1_secondHalfEdge[2];
}
else // we look from smaller cell = unique interface
{
// proceed with numbering according to Aavatsmark, seen from cell i
irregularInterfaceMap_[intersectionID].T1_[first] = T1;
irregularInterfaceMap_[intersectionID].T1_[second] = T11_secondHalfEdge;
irregularInterfaceMap_[intersectionID].T1_[second] = T1_secondHalfEdge;
}
irregularInterfaceMap_[intersectionID].globalPos3_[0] = globalPos3;
......@@ -228,10 +259,24 @@ public:
return;
}
//! Provides acess to stored Mpfa Data on an intersection for both half-edges
/** The method gets information of both interaction regions (Transmissitivity
* as well as details of the 3rd cells of both regions) into a storage container.
* The key to each element is the index of the intersection, seen from the smaller
* cell (only this is unique).
* If we arrive from the "wrong" (i.e. non-unique) direction, we invert fluxes.
*
* \param irregularIs The current irregular intersection
* \param secondHalfEdgeIntersectionIt Iterator to the intersection connecting the second interaction region
* \param T1 Transmissitivity matrix for flux calculations: unique interaction region
* \param T1_secondHalfEdge Second transmissitivity matrix for flux calculations for non-unique interaction region
* \param globalPos3 The position of the 3rd cell of the first interaction region
* \param globalIdx3 The index of the 3rd cell of the first interaction region
*/
int getMpfaData(const Intersection& irregularIs,
IntersectionIterator& secondHalfEdgeIntersectionIt,
T1ransmissivityMatrix& T1,
T1ransmissivityMatrix& T1_secondHalfEdge,
TransmissivityMatrix& T1,
TransmissivityMatrix& T1_secondHalfEdge,
GlobalPosition& globalPos3,
int& globalIdx3)
{
......
......@@ -19,7 +19,7 @@
/*!
* \file
*
* \brief test problem for the sequential 2p2c model
* \brief Test problem for the adaptive sequential 2p2c model
*/
#ifndef DUMUX_TEST_ADAPTIVE_2P2C_PROBLEM_HH
#define DUMUX_TEST_ADAPTIVE_2P2C_PROBLEM_HH
......@@ -102,6 +102,7 @@ SET_INT_PROP(Adaptive2p2c, PressureFormulation,
/*!
* \ingroup Adaptive2p2cs
* \ingroup IMPETtests
*/
template<class TypeTag = TTAG(Adaptive2p2c)>