From 72076121ac329a08fd926571304f2dda206e4163 Mon Sep 17 00:00:00 2001
From: Markus Wolff <markus.wolff@twt-gmbh.de>
Date: Fri, 17 Feb 2012 15:36:08 +0000
Subject: [PATCH] improved doxygen documentation of 2p standard finite volume
 pressure model

git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@7813 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
 .../decoupled/2p/diffusion/fv/fvpressure2p.hh | 212 +++++++++++++-----
 .../decoupled/2p/diffusion/fv/fvvelocity2p.hh |  65 ++++--
 2 files changed, 203 insertions(+), 74 deletions(-)

diff --git a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
index 33ff2a1c71..67f6371d62 100644
--- a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
@@ -30,27 +30,60 @@
 
 /**
  * @file
- * @brief  Finite Volume Discretization of a pressure equation.
- * @author Bernd Flemisch, Jochen Fritz, Markus Wolff
+ * @brief  Finite Volume Discretization of a two-phase flow pressure equation.
+ * @author Markus Wolff, Bernd Flemisch, Jochen Fritz,
  */
 
 namespace Dumux
 {
-//! \ingroup FV2p
-//! \brief Finite Volume discretization of the pressure equation of the sequential IMPES Model.
-/*! Provides a Finite Volume implementation for the evaluation
- * of equations of the form
- * \f[\text{div}\, \boldsymbol{v}_{total} = q.\f]
- * The definition of the total velocity \f$\boldsymbol{v}_total\f$ depends on the kind of pressure chosen. This could be a wetting (w) phase pressure leading to
- * \f[ - \text{div}\,  \left[\lambda \boldsymbol{K} \left(\text{grad}\, p_w + f_n \text{grad}\, p_c + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q, \f]
- * a non-wetting (n) phase pressure yielding
- * \f[ - \text{div}\,  \left[\lambda \boldsymbol{K}  \left(\text{grad}\, p_n - f_w \text{grad}\, p_c + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q, \f]
- * or a global pressure leading to
- * \f[ - \text{div}\, \left[\lambda \boldsymbol{K} \left(\text{grad}\, p_{global} + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q.\f]
- *  Here, \f$p\f$ denotes a pressure, \f$\boldsymbol{K}\f$ the absolute permeability, \f$\lambda\f$ the total mobility, possibly depending on the
- * saturation,\f$f\f$ the fractional flow function of a phase, \f$\rho\f$ a phase density, \f$g\f$ the gravity constant and \f$q\f$ the source term.
- * For all cases, \f$p = p_D\f$ on \f$\Gamma_{Neumann}\f$, and \f$\boldsymbol{v}_{total}  = q_N\f$
- * on \f$\Gamma_{Dirichlet}\f$.
+//! \ingroup FVPressure2p
+/*!  \brief Finite Volume discretization of a two-phase flow pressure equation of the sequential IMPES Model.
+ *
+ * This class provides a finite volume (FV) implementation for solving equations of the form
+ * \f[
+ * \phi \left( \rho_w  \frac{\partial S_w}{\partial t} + \rho_n \frac{\partial S_n}{\partial t}\right) + \text{div}\, \boldsymbol{v}_{total} = q.
+ * \f]
+ * The definition of the total velocity \f$\boldsymbol{v}_{total}\f$ depends on the choice of the primary pressure variable.
+ * Further, fluids can be assumed to be compressible or incompressible (Property: <tt>EnableCompressibility</tt>).
+ * In the incompressible case a wetting (\f$ w \f$) phase pressure as primary variable leads to
+ *
+ * \f[
+ * - \text{div}\,  \left[\lambda \boldsymbol K \left(\text{grad}\, p_w + f_n \text{grad}\, p_c + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q,
+ * \f]
+ *
+ * a non-wetting (\f$ n \f$) phase pressure yields
+ * \f[
+ *  - \text{div}\,  \left[\lambda \boldsymbol K  \left(\text{grad}\, p_n - f_w \text{grad}\, p_c + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q,
+ *  \f]
+ * and a global pressure leads to
+ * \f[
+ * - \text{div}\, \left[\lambda \boldsymbol K \left(\text{grad}\, p_{global} + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q.
+ * \f]
+ * Here, \f$ p_\alpha \f$ is a phase pressure, \f$ p_ {global} \f$ the global pressure of a classical fractional flow formulation
+ * (see e.g. P. Binning and M. A. Celia, “Practical implementation of the fractional flow approach to multi-phase flow simulation,” Advances in water resources, vol. 22, no. 5, pp. 461-478, 1999.),
+ * \f$ p_c = p_n - p_w \f$ is the capillary pressure, \f$ \boldsymbol K \f$ the absolute permeability, \f$ \lambda = \lambda_w +  \lambda_n \f$ the total mobility depending on the
+ * saturation (\f$ \lambda_\alpha = k_{r_\alpha} / \mu_\alpha \f$),\f$ f_\alpha = \lambda_\alpha / \lambda \f$ the fractional flow function of a phase,
+ * \f$ \rho_\alpha \f$ a phase density, \f$ g \f$ the gravity constant and \f$ q \f$ the source term.
+ *
+ * For all cases, \f$ p = p_D \f$ on \f$ \Gamma_{Neumann} \f$, and \f$ \boldsymbol{v}_{total}  = q_N \f$
+ * on \f$ \Gamma_{Dirichlet} \f$.
+ *
+ * The slightly compressible case is only implemented for phase pressures! In this case for a wetting (\f$ w \f$) phase pressure as primary variable the equations are formulated as
+ * \f[
+ * \phi \left( \rho_w  \frac{\partial S_w}{\partial t} + \rho_n \frac{\partial S_n}{\partial t}\right) - \text{div}\,  \left[\lambda \boldsymbol{K} \left(\text{grad}\, p_w + f_n \text{grad}\, p_c + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q,
+ * \f]
+ * and for a non-wetting (\f$ n \f$) phase pressure as
+ *  \f[
+ *  \phi \left( \rho_w  \frac{\partial S_w}{\partial t} + \rho_n \frac{\partial S_n}{\partial t}\right) - \text{div}\,  \left[\lambda \boldsymbol{K}  \left(\text{grad}\, p_n - f_w \text{grad}\, p_c + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right)\right] = q,
+ *  \f]
+ * In this slightly compressible case the following definitions are valid:  \f$ \lambda = \rho_w \lambda_w + \rho_n \lambda_n \f$, \f$ f_\alpha = (\rho_\alpha \lambda_\alpha) / \lambda \f$
+ * This model assumes that temporal changes in density are very small and thus terms of temporal derivatives are negligible in the pressure equation.
+ * Depending on the formulation the terms including time derivatives of saturations are simplified by inserting  \f$ S_w + S_n = 1 \f$.
+ *
+ *  In the IMPES models the default setting is:
+ *
+ *      - formulation: \f$ p_w-S_w \f$ (Property: <tt>Formulation</tt> defined as <tt>DecoupledTwoPCommonIndices::pwSw</tt>)
+ *      - compressibility: disabled (Property: <tt>EnableCompressibility</tt> set to <tt>false</tt>)
  *
  * \tparam TypeTag The Type Tag
  */
@@ -106,26 +139,41 @@ template<class TypeTag> class FVPressure2P: public FVPressure<TypeTag>
     typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
     typedef Dune::FieldMatrix<Scalar, dim, dim> FieldMatrix;
 
-public:
-
+protected:
+    //! \cond \private
     enum
     {
         rhs = ParentType::rhs, matrix = ParentType::matrix
     };
+    //! \endcond
+
+public:
+
+    //! \cond \private
 
-    void getSource(Dune::FieldVector<Scalar, 2>&, const Element&, const CellData&, const bool);
+    // Function which calculates the source entry
+    void getSource(Dune::FieldVector<Scalar, 2>& entry, const Element& element, const CellData& cellData, const bool first);
 
-    void getStorage(Dune::FieldVector<Scalar, 2>&, const Element&, const CellData&, const bool);
+    // Function which calculates the storage entry
+    void getStorage(Dune::FieldVector<Scalar, 2>& entry, const Element& element, const CellData& cellData, const bool first);
 
-    void getFlux(Dune::FieldVector<Scalar, 2>&, const Intersection&, const CellData&, const bool);
+    // Function which calculates the flux entry
+    void getFlux(Dune::FieldVector<Scalar, 2>& entry, const Intersection& intersection, const CellData& cellData, const bool first);
 
-    void getFluxOnBoundary(Dune::FieldVector<Scalar, 2>&,
-    const Intersection&, const CellData&, const bool);
+    // Function which calculates the boundary flux entry
+    void getFluxOnBoundary(Dune::FieldVector<Scalar, 2>& entry,
+    const Intersection& intersection, const CellData& cellData, const bool first);
 
-    //! updates and stores constitutive relations
+    // updates and stores constitutive relations
     void updateMaterialLaws();
+    //! \endcond
 
-    //! \copydoc Dumux::FVPressure1P::initialize()
+    /*! \brief Initializes the pressure model
+     *
+     * \copydoc ParentType::initialize()
+     *
+     * \param solveTwice indicates if more than one iteration is allowed to get an initial pressure solution
+     */
     void initialize(bool solveTwice = true)
     {
         ParentType::initialize();
@@ -169,13 +217,17 @@ public:
         return;
     }
 
-    //! updates the pressure field (analog to update function in Dumux::IMPET)
+    /*! \brief Pressure update
+     *
+     * \copydoc ParentType::update()
+     *
+     */
     void update()
     {
+        timeStep_ = problem_.timeManager().timeStepSize();
         //error bounds for error term for incompressible models to correct unphysical saturation over/undershoots due to saturation transport
         if (!compressibility_)
         {
-            timeStep_ = problem_.timeManager().timeStepSize();
             maxError_ = 0.0;
             int size = problem_.gridView().size(0);
             for (int i = 0; i < size; i++)
@@ -208,6 +260,9 @@ public:
         return;
     }
 
+    /*! \brief Globally stores the pressure solution
+     *
+     */
     void storePressureSolution()
     {
         int size = problem_.gridView().size(0);
@@ -219,6 +274,13 @@ public:
         }
     }
 
+    /*! \brief Stores the pressure solution of a cell
+     *
+     * Calculates secondary pressure variables and stores pressures.
+     *
+     * \param globalIdx Global cell index
+     * \param cellData A CellData object
+     */
     void storePressureSolution(int globalIdx, CellData& cellData)
     {
         switch (pressureType_)
@@ -251,7 +313,15 @@ public:
         }
     }
 
-    //! \copydoc Dumux::FVPressure1P::addOutputVtkFields(MultiWriter &writer)
+    /*! \brief Adds pressure output to the output file
+     *
+     * Adds the phase pressures or a global pressure (depending on the formulation) as well as the capillary pressure to the output.
+     * In the compressible case additionally density and viscosity are added.
+     *
+     * \tparam MultiWriter Class defining the output writer
+     * \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
+     *
+     */
     template<class MultiWriter>
     void addOutputVtkFields(MultiWriter &writer)
     {
@@ -368,9 +438,14 @@ private:
     static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation); //!< gives kind of saturation used (\f$S_w\f$, \f$S_n\f$)
 };
 
-//!function which calculates the source entry
+/*! \brief Function which calculates the source entry
+ *
+ * \copydoc FVPressure::getSource(Dune::FieldVector<Scalar,2>&,const Element&,const CellData&,const bool)
+ *
+ * Source of each fluid phase has to be added as mass flux (\f$\text{kg}/(\text{m}^3 \text{s}\f$).
+ */
 template<class TypeTag>
-void FVPressure2P<TypeTag>::getSource(Dune::FieldVector<Scalar, 2>& entries, const Element& element
+void FVPressure2P<TypeTag>::getSource(Dune::FieldVector<Scalar, 2>& entry, const Element& element
         , const CellData& cellData, const bool first)
 {
     // cell volume, assume linear map here
@@ -386,14 +461,34 @@ void FVPressure2P<TypeTag>::getSource(Dune::FieldVector<Scalar, 2>& entries, con
         sourcePhase[nPhaseIdx] /= density_[nPhaseIdx];
     }
 
-    entries[rhs] = volume * (sourcePhase[wPhaseIdx] + sourcePhase[nPhaseIdx]);
+    entry[rhs] = volume * (sourcePhase[wPhaseIdx] + sourcePhase[nPhaseIdx]);
 
     return;
 }
 
-//!function which calculates the storage entry
+/** \brief Function which calculates the storage entry
+ *
+ * \copydoc FVPressure::getStorage(Dune::FieldVector<Scalar,2>&,const Element&,const CellData&,const bool)
+ *
+ * If compressibility is enabled this functions calculates the term
+ * \f[
+ *      \phi \sum_\alpha \rho_\alpha \frac{\partial S_\alpha}{\partial t} V
+ * \f]
+ *
+ * In the incompressible case an volume correction term is calculated which corrects for unphysical saturation overshoots/undershoots.
+ * These can occur if the estimated time step for the explicit transport was too large. Correction by an artificial source term allows to correct
+ * this errors due to wrong time-stepping without losing mass conservation. The error term looks as follows:
+ * \f[
+ *  q_{error} = \begin{cases}
+ *          S < 0 & a_{error} \frac{S}{\Delta t} V \\
+ *          S > 1 & a_{error} \frac{(S - 1)}{\Delta t} V \\
+ *          0 \le S \le 1 & 0
+ *      \end{cases}
+ *  \f]
+ *  where \f$a_{error}\f$ is a weighting factor (default: \f$a_{error} = 0.5\f$)
+ */
 template<class TypeTag>
-void FVPressure2P<TypeTag>::getStorage(Dune::FieldVector<Scalar, 2>& entries, const Element& element
+void FVPressure2P<TypeTag>::getStorage(Dune::FieldVector<Scalar, 2>& entry, const Element& element
         , const CellData& cellData, const bool first)
 {
     //volume correction due to density differences
@@ -408,13 +503,13 @@ void FVPressure2P<TypeTag>::getStorage(Dune::FieldVector<Scalar, 2>& entries, co
         {
         case Sw:
         {
-            entries[rhs] = -(cellData.volumeCorrection() * porosity * volume
+            entry[rhs] = -(cellData.volumeCorrection()/timeStep_ * porosity * volume
                     * (cellData.density(wPhaseIdx) - cellData.density(nPhaseIdx)));
             break;
         }
         case Sn:
         {
-            entries[rhs] = -(cellData.volumeCorrection() * porosity * volume
+            entry[rhs] = -(cellData.volumeCorrection()/timeStep_ * porosity * volume
                     * (cellData.density(nPhaseIdx) - cellData.density(wPhaseIdx)));
             break;
         }
@@ -445,16 +540,20 @@ void FVPressure2P<TypeTag>::getStorage(Dune::FieldVector<Scalar, 2>& entries, co
 
         if ((errorAbs*timeStep_ > 1e-6) && (errorAbs > ErrorTermLowerBound_ * maxError_) && (!problem_.timeManager().willBeFinished()))
         {
-            entries[rhs] = ErrorTermFactor_ * error * volume;
+            entry[rhs] = ErrorTermFactor_ * error * volume;
         }
     }
 
     return;
 }
 
-//!function which calculates internal flux entries
+/*! \brief Function which calculates the flux entry
+ *
+ * \copydoc FVPressure::getFlux(Dune::FieldVector<Scalar,2>&,const Intersection&,const CellData&,const bool)
+ *
+ */
 template<class TypeTag>
-void FVPressure2P<TypeTag>::getFlux(Dune::FieldVector<Scalar, 2>& entries, const Intersection& intersection
+void FVPressure2P<TypeTag>::getFlux(Dune::FieldVector<Scalar, 2>& entry, const Intersection& intersection
         , const CellData& cellDataI, const bool first)
 {
     ElementPointer elementI = intersection.inside();
@@ -562,10 +661,10 @@ void FVPressure2P<TypeTag>::getFlux(Dune::FieldVector<Scalar, 2>& entries, const
     }
 
     //calculate current matrix entry
-    entries[matrix] = (lambdaW + lambdaNW) * ((permeability * unitOuterNormal) / dist) * faceArea;
+    entry[matrix] = (lambdaW + lambdaNW) * ((permeability * unitOuterNormal) / dist) * faceArea;
 
     //calculate right hand side
-    entries[rhs] = (lambdaW * density_[wPhaseIdx] + lambdaNW * density_[nPhaseIdx]) * (permeability * gravity_) * faceArea;
+    entry[rhs] = (lambdaW * density_[wPhaseIdx] + lambdaNW * density_[nPhaseIdx]) * (permeability * gravity_) * faceArea;
 
     switch (pressureType_)
     {
@@ -576,7 +675,7 @@ void FVPressure2P<TypeTag>::getFlux(Dune::FieldVector<Scalar, 2>& entries, const
         pCGradient *= (pcI - pcJ) / dist;
 
         //add capillary pressure term to right hand side
-        entries[rhs] += 0.5 * (lambdaNWI + lambdaNWJ) * (permeability * pCGradient) * faceArea;
+        entry[rhs] += 0.5 * (lambdaNWI + lambdaNWJ) * (permeability * pCGradient) * faceArea;
         break;
     }
     case pn:
@@ -586,7 +685,7 @@ void FVPressure2P<TypeTag>::getFlux(Dune::FieldVector<Scalar, 2>& entries, const
         pCGradient *= (pcI - pcJ) / dist;
 
         //add capillary pressure term to right hand side
-        entries[rhs] -= 0.5 * (lambdaWI + lambdaWJ) * (permeability * pCGradient) * faceArea;
+        entry[rhs] -= 0.5 * (lambdaWI + lambdaWJ) * (permeability * pCGradient) * faceArea;
         break;
     }
     }
@@ -594,9 +693,15 @@ void FVPressure2P<TypeTag>::getFlux(Dune::FieldVector<Scalar, 2>& entries, const
     return;
 }
 
-//!function which calculates internal flux entries
+/*! \brief Function which calculates the flux entry at a boundary
+ *
+ * \copydoc FVPressure::getFluxOnBoundary(Dune::FieldVector<Scalar,2>&,const Intersection&,const CellData&,const bool)
+ *
+ * Dirichlet boundary condition is a pressure depending on the formulation (\f$p_w\f$ (default), \f$p_n\f$, \f$p_{global}\f$),
+ * Neumann boundary condition are the phase mass fluxes (\f$q_w\f$ and \f$q_n\f$, [\f$\text{kg}/(\text{m}^2 \text{s}\f$])
+ */
 template<class TypeTag>
-void FVPressure2P<TypeTag>::getFluxOnBoundary(Dune::FieldVector<Scalar, 2>& entries,
+void FVPressure2P<TypeTag>::getFluxOnBoundary(Dune::FieldVector<Scalar, 2>& entry,
 const Intersection& intersection, const CellData& cellData, const bool first)
 {
     ElementPointer element = intersection.inside();
@@ -797,11 +902,11 @@ const Intersection& intersection, const CellData& cellData, const bool first)
         }
 
         //calculate current matrix entry
-        entries[matrix] = (lambdaW + lambdaNW) * ((permeability * unitOuterNormal) / dist) * faceArea;
-        entries[rhs] = entries[matrix] * pressBound;
+        entry[matrix] = (lambdaW + lambdaNW) * ((permeability * unitOuterNormal) / dist) * faceArea;
+        entry[rhs] = entry[matrix] * pressBound;
 
         //calculate right hand side
-        entries[rhs] -= (lambdaW * density_[wPhaseIdx] + lambdaNW * density_[nPhaseIdx]) * (permeability * gravity_)
+        entry[rhs] -= (lambdaW * density_[wPhaseIdx] + lambdaNW * density_[nPhaseIdx]) * (permeability * gravity_)
                 * faceArea;
 
         switch (pressureType_)
@@ -813,7 +918,7 @@ const Intersection& intersection, const CellData& cellData, const bool first)
             pCGradient *= (pcI - pcBound) / dist;
 
             //add capillary pressure term to right hand side
-            entries[rhs] -= 0.5 * (lambdaNWI + lambdaNWBound) * (permeability * pCGradient) * faceArea;
+            entry[rhs] -= 0.5 * (lambdaNWI + lambdaNWBound) * (permeability * pCGradient) * faceArea;
             break;
         }
         case pn:
@@ -823,7 +928,7 @@ const Intersection& intersection, const CellData& cellData, const bool first)
             pCGradient *= (pcI - pcBound) / dist;
 
             //add capillary pressure term to right hand side
-            entries[rhs] += 0.5 * (lambdaWI + lambdaWBound) * (permeability * pCGradient) * faceArea;
+            entry[rhs] += 0.5 * (lambdaWI + lambdaWBound) * (permeability * pCGradient) * faceArea;
             break;
         }
         }
@@ -838,7 +943,7 @@ const Intersection& intersection, const CellData& cellData, const bool first)
             boundValues[wPhaseIdx] /= density_[wPhaseIdx];
             boundValues[nPhaseIdx] /= density_[nPhaseIdx];
         }
-        entries[rhs] = -(boundValues[wPhaseIdx] + boundValues[nPhaseIdx]) * faceArea;
+        entry[rhs] = -(boundValues[wPhaseIdx] + boundValues[nPhaseIdx]) * faceArea;
     }
     else
     {
@@ -848,7 +953,10 @@ const Intersection& intersection, const CellData& cellData, const bool first)
     return;
 }
 
-//!constitutive functions are updated once if new saturations are calculated and stored in the variables object
+/*! \brief Updates constitutive relations and stores them in the variable class
+ *
+ * Stores mobility, fractional flow function and capillary pressure for all grid cells. In the compressible case additionally the densities and viscosities are stored.
+ */
 template<class TypeTag>
 void FVPressure2P<TypeTag>::updateMaterialLaws()
 {
diff --git a/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh b/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
index 58efa62322..83c2e20ffb 100644
--- a/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
+++ b/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
@@ -33,19 +33,23 @@
 
 namespace Dumux
 {
-//! \ingroup FV2p
-//! \brief Determines the velocity from a finite volume solution of the  pressure equation of the sequential Model (IMPES).
-/*! Calculates phase velocities or total velocity from a known pressure field in context of a Finite Volume implementation for the evaluation
- * of equations of the form
- * \f[\text{div}\, \boldsymbol{v}_{total} = q.\f]
+//! \ingroup FVPressure2p
+//! \brief Determines the velocity from a finite volume solution of the  pressure equation of a sequential model (IMPES).
+/*! Calculates phase velocities or total velocity from a known pressure field applying a finite volume discretization.
  * The wetting or the non-wetting phase pressure, or the global pressure has to be given as piecewise constant cell values.
  * The phase velocities are calculated following  Darcy's law as
- * \f[\boldsymbol{v}_\alpha = \lambda_\alpha \boldsymbol{K} \left(\text{grad}\, p_\alpha + \rho_\alpha g  \text{grad}\, z\right),\f]
- * where \f$p_\alpha\f$ denotes the pressure of phase \f$_\alpha\f$ (wetting or non-wetting), \f$\boldsymbol{K}\f$ the absolute permeability, \f$\lambda_\alpha\f$ the phase mobility, \f$\rho_\alpha\f$ the phase density and \f$g\f$ the gravity constant.
+ * \f[
+ * \boldsymbol v_\alpha = \lambda_\alpha \boldsymbol K \left(\text{grad}\, p_\alpha + \rho_\alpha g  \text{grad}\, z \right),
+ * \f]
+ * where \f$ p_\alpha \f$ denotes the pressure of phase \f$ _\alpha \f$ (wetting or non-wetting), \f$ \boldsymbol K \f$ the absolute permeability, \f$ \lambda_\alpha \f$ the phase mobility, \f$ \rho_\alpha \f$ the phase density and \f$ g \f$ the gravity constant.
  * The total velocity is either calculated as sum of the phase velocities
- * \f[\boldsymbol{v}_{total} = \boldsymbol{v}_{wetting}+\boldsymbol{v}_{non-wetting},\f]
+ * \f[
+ * \boldsymbol v_{total} = \boldsymbol v_{wetting}+\boldsymbol v_{non-wetting},
+ * \f]
  * or with a given global pressure
- * \f[\boldsymbol{v}_{total} = \lambda_{total} \boldsymbol{K} \left(\text{grad}\, p_{global} + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right).\f]
+ * \f[
+ * \boldsymbol v_{total} = \lambda_{total} \boldsymbol K \left(\text{grad}\, p_{global} + \sum f_\alpha \rho_\alpha g  \text{grad}\, z\right).
+ * \f]
  *
  * \tparam TypeTag The Type Tag
  */
@@ -107,9 +111,8 @@ class FVVelocity2P
     typedef Dune::FieldMatrix<Scalar, dim, dim> FieldMatrix;
 
 public:
-    //! Constructs a FVVelocity2P object
-    /*!
-     * \param problem a problem class object
+    /*! \brief Constructs a FVVelocity2P object
+     * \param problem A Problem class object
      */
     FVVelocity2P(Problem& problem) :
             problem_(problem), gravity_(problem.gravity())
@@ -140,25 +143,29 @@ public:
         }
     }
 
-    //! Calculate the velocity.
-    /*!
-     *
-     *  Given the piecewise constant pressure \f$p\f$,
-     *  this method calculates the velocity
-     *  The method is needed in the IMPES (Implicit Pressure Explicit Saturation) algorithm which is used for a fractional flow formulation
-     *  to provide the velocity field required for the solution of the saturation equation.
-     */
+    //! \cond \private
     void calculateVelocity(const Intersection&, CellData&);
 
     void calculateVelocityOnBoundary(const Intersection&, CellData&);
+    //! \endcond
 
+    /*! \brief Indicates if velocity is reconstructed in the pressure step or in the transport step
+     *
+     * Returns true (In the standard finite volume discretization the velocity is calculated during the saturation transport.)
+     */
     bool calculateVelocityInTransport()
     {
         return true;
     }
 
-    //! \brief Write data files
-    /*  \param name file name */
+    /*! \brief Adds velocity output to the output file
+     *
+     * Adds the phase velocities or a total velocity (depending on the formulation) to the output.
+     *
+     * \tparam MultiWriter Class defining the output writer
+     * \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
+     *
+     */
     template<class MultiWriter>
     void addOutputVtkFields(MultiWriter &writer)
     {
@@ -264,6 +271,13 @@ private:
     static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation); //!< gives kind of saturation used (\f$S_w\f$, \f$S_n\f$)
 };
 
+/*! \brief Calculates the velocity at a cell-cell interface.
+*
+* Calculates the velocity at a cell-cell interface from a given pressure field.
+*
+* \param intersection Intersection of two grid cells
+* \param CellData Object containing all model relevant cell data
+*/
 template<class TypeTag>
 void FVVelocity2P<TypeTag>::calculateVelocity(const Intersection& intersection, CellData& cellDataI)
 {
@@ -429,6 +443,13 @@ void FVVelocity2P<TypeTag>::calculateVelocity(const Intersection& intersection,
     return;
 }
 
+/*! \brief Calculates the velocity at a boundary.
+*
+* Calculates the velocity at a boundary from a given pressure field.
+*
+* \param intersection Boundary intersection
+* \param CellData Object containing all model relevant cell data
+*/
 template<class TypeTag>
 void FVVelocity2P<TypeTag>::calculateVelocityOnBoundary(const Intersection& intersection, CellData& cellDataI)
 {
-- 
GitLab