diff --git a/doc/handbook/ModelDescriptions/1p2cboxmodel.tex b/doc/handbook/ModelDescriptions/1p2cboxmodel.tex
index adf7432c7487edfe11516146989d648b14953b74..37e3243117d42d4265219b03ba0f8fd96a2c9f4f 100644
--- a/doc/handbook/ModelDescriptions/1p2cboxmodel.tex
+++ b/doc/handbook/ModelDescriptions/1p2cboxmodel.tex
@@ -4,14 +4,13 @@
% file instead!! %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\-This model implements a one-\/phase flow of a compressible fluid, that consists of two components, using a standard \-Darcy approach as the equation for the conservation of momentum\-: \[ v_{D} = - \frac{\textbf K}{\mu} \left(\text{grad} p - \varrho {\textbf g} \right) \]
-Adaption of the BOX scheme to the one-\/phase two-\/component flow model. This model implements an one-\/phase flow of an incompressible fluid, that consists of two components, using a standard Darcy approach as the equation for the conservation of momentum: \[ v_{D} = - \frac{K}{\mu} \left(\text{grad} p - \varrho g \right) \]
+\-Gravity can be enabled or disabled via the property system. \-By inserting this into the continuity equation, one gets \[ \Phi \frac{\partial \varrho}{\partial t} - \text{div} \left\{ \varrho \frac{\textbf K}{\mu} \left(\text{grad}\, p - \varrho {\textbf g} \right) \right\} = q \;, \]
-Gravity can be enabled or disabled via the Property system. By inserting this into the continuity equation, one gets \[ - \text{div} \left\{ \varrho \frac{K}{\mu} \left(\text{grad} p - \varrho g \right) \right\} = q \;, \]
+\-The transport of the components is described by the following equation\-: \[ \Phi \frac{ \partial \varrho x}{\partial t} - \text{div} \left( \varrho \frac{{\textbf K} x}{\mu} \left( \text{grad}\, p - \varrho {\textbf g} \right) + \varrho \tau \Phi D \text{grad} x \right) = q. \]
-The transport of the components is described by the following equation: \[ \Phi \varrho \frac{ \partial x}{\partial t} - \text{div} \left( \varrho \frac{K x}{\mu} \left( \text{grad} p - \varrho g \right) + \varrho \tau \Phi D \text{grad} x \right) = q. \]
+\-All equations are discretized using a fully-\/coupled vertex-\/centered finite volume (box) scheme as spatial and the implicit \-Euler method as time discretization.
-All equations are discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit Euler method as time discretization.
-
-The primary variables are the pressure $p$ and the mole fraction of dissolved component $x$.
+\-The primary variables are the pressure $p$ and the mole or mass fraction of dissolved component $x$.
diff --git a/doc/handbook/ModelDescriptions/1pboxmodel.tex b/doc/handbook/ModelDescriptions/1pboxmodel.tex
index 91b93512e6dd27fad827c6c03a1075f78ad1d3ca..51bf610363abd045acd527fa2fcdeeb48c8c4bdf 100644
--- a/doc/handbook/ModelDescriptions/1pboxmodel.tex
+++ b/doc/handbook/ModelDescriptions/1pboxmodel.tex
@@ -4,6 +4,5 @@
% file instead!! %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-Adaption of the BOX scheme to the single phase isothermal flow model. Single phase compressible isothermal flow model, \begin{align*} \phi \frac{\partial \varrho}{\partial t} + \vec{\nabla} \cdot (- \varrho \frac{\bar{\bar{K}}}{\mu} ( \nabla p -\varrho \vec{g})) = q, \end{align*} discretized using a vertex centered finite volume (box) scheme as spatial and the implicit Euler method as time discretization. Of course, the model can also be used for incompressible single phase flow modeling, if in the problem file a fluid with constant density is chosen.
+\-Single-\/phase compressible isothermal flow model, \begin{align*} \phi \frac{\partial \varrho}{\partial t} + \text{div} (- \varrho \frac{\textbf K}{\mu} ( \text{grad}\, p -\varrho {\textbf g})) = q, \end{align*} discretized using a vertex-\/centered finite volume (box) scheme as spatial and the implicit \-Euler method as time discretization. \-Of course, the model can also be used for incompressible single phase flow modeling, if a fluid with constant density is chosen in the problem file.
diff --git a/doc/handbook/ModelDescriptions/2p2cboxmodel.tex b/doc/handbook/ModelDescriptions/2p2cboxmodel.tex
index fb378aa01b284a861a17596e24698a51abd71f8b..6e39c9d747e45209b70b8a1d2eebe8c7f7ab79d7 100644
--- a/doc/handbook/ModelDescriptions/2p2cboxmodel.tex
+++ b/doc/handbook/ModelDescriptions/2p2cboxmodel.tex
@@ -4,17 +4,16 @@
% file instead!! %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\-This model implements two-\/phase two-\/component flow of two compressible and partially miscible fluids $\alpha \in \{ w, n \}$ composed of the two components $\kappa \in \{ w, a \}$. \-The standard multiphase \-Darcy approach is used as the equation for the conservation of momentum\-: \[ v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left(\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g} \right) \]
-Adaption of the BOX scheme to the two-\/phase two-\/component flow model. This model implements two-\/phase two-\/component flow of two compressible and partially miscible fluids $\alpha \in \{ w, n \}$ composed of the two components $\kappa \in \{ w, a \}$. The standard multiphase Darcy approach is used as the equation for the conservation of momentum: \[ v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left(\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g} \right) \]
+\-By inserting this into the equations for the conservation of the components, one gets one transport equation for each component \begin{eqnarray} && \phi \frac{\partial (\sum_\alpha \varrho_\alpha X_\alpha^\kappa S_\alpha )} {\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha X_\alpha^\kappa \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} (\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g}) \right\} \nonumber \\ \nonumber \\ &-& \sum_\alpha \text{div} \left\{{\bf D}_{\alpha, pm}^\kappa \varrho_{\alpha} \text{grad}\, X^\kappa_{\alpha} \right\} - \sum_\alpha q_\alpha^\kappa = 0 \qquad \kappa \in \{w, a\} \, , \alpha \in \{w, g\} \end{eqnarray}
-By inserting this into the equations for the conservation of the components, one gets one transport equation for each component \begin{eqnarray} && \phi \frac{\partial (\sum_\alpha \varrho_\alpha X_\alpha^\kappa S_\alpha )} {\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha X_\alpha^\kappa \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} (\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g}) \right\} \nonumber \\ \nonumber \\ &-& \sum_\alpha \text{div} \left\{{\bf D_{\alpha, pm}^\kappa} \varrho_{\alpha} \text{grad}\, X^\kappa_{\alpha} \right\} - \sum_\alpha q_\alpha^\kappa = 0 \qquad \kappa \in \{w, a\} \, , \alpha \in \{w, g\} \end{eqnarray}
+\-This is discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit \-Euler method as temporal discretization.
-This is discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit Euler method as temporal discretization.
-
-By using constitutive relations for the capillary pressure $p_c = p_n - p_w$ and relative permeability $k_{r\alpha}$ and taking advantage of the fact that $S_w + S_n = 1$ and $X^\kappa_w + X^\kappa_n = 1$, the number of unknowns can be reduced to two. The used primary variables are, like in the two-\/phase model, either $p_w$ and $S_n$ or $p_n$ and $S_w$. The formulation which ought to be used can be specified by setting the {\ttfamily Formulation} property to either TwoPTwoCIndices::pWsN or TwoPTwoCIndices::pNsW. By default, the model uses $p_w$ and $S_n$. Moreover, the second primary variable depends on the phase state, since a primary variable switch is included. The phase state is stored for all nodes of the system. Following cases can be distinguished:
+\-By using constitutive relations for the capillary pressure $p_c = p_n - p_w$ and relative permeability $k_{r\alpha}$ and taking advantage of the fact that $S_w + S_n = 1$ and $X^\kappa_w + X^\kappa_n = 1$, the number of unknowns can be reduced to two. \-The used primary variables are, like in the two-\/phase model, either $p_w$ and $S_n$ or $p_n$ and $S_w$. \-The formulation which ought to be used can be specified by setting the {\ttfamily \-Formulation} property to either \-Two\-P\-Two\-C\-Indices\-::p\-Ws\-N or \-Two\-P\-Two\-C\-Indices\-::p\-Ns\-W. \-By default, the model uses $p_w$ and $S_n$. \-Moreover, the second primary variable depends on the phase state, since a primary variable switch is included. \-The phase state is stored for all nodes of the system. \-Following cases can be distinguished\-:
\begin{itemize}
-\item Both phases are present: The saturation is used (either $S_n$ or $S_w$, dependent on the chosen {\ttfamily Formulation}), as long as $ 0 < S_\alpha < 1$.
-\item Only wetting phase is present: The mass fraction of, e.g., air in the wetting phase $X^a_w$ is used, as long as the maximum mass fraction is not exceeded ( $X^a_w<X^a_{w,max}$)
-\item Only non-\/wetting phase is present: The mass fraction of, e.g., water in the non-\/wetting phase, $X^w_n$, is used, as long as the maximum mass fraction is not exceeded ( $X^w_n<X^w_{n,max}$)
+\item \-Both phases are present\-: \-The saturation is used (either $S_n$ or $S_w$, dependent on the chosen {\ttfamily \-Formulation}), as long as $ 0 < S_\alpha < 1$.
+\item \-Only wetting phase is present\-: \-The mass fraction of, e.\-g., air in the wetting phase $X^a_w$ is used, as long as the maximum mass fraction is not exceeded ( $X^a_w<X^a_{w,max}$)
+\item \-Only non-\/wetting phase is present\-: \-The mass fraction of, e.\-g., water in the non-\/wetting phase, $X^w_n$, is used, as long as the maximum mass fraction is not exceeded ( $X^w_n<X^w_{n,max}$)
\end{itemize}
diff --git a/doc/handbook/ModelDescriptions/2p2cniboxmodel.tex b/doc/handbook/ModelDescriptions/2p2cniboxmodel.tex
index 1116e4512c38616b355f58e66b75f8ab75dbbc7d..c77d992b8a16b6e6d8a52b5aea6feb401af64d56 100644
--- a/doc/handbook/ModelDescriptions/2p2cniboxmodel.tex
+++ b/doc/handbook/ModelDescriptions/2p2cniboxmodel.tex
@@ -4,15 +4,14 @@
% file instead!! %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\-This model implements a non-\/isothermal two-\/phase flow of two compressible and partly miscible fluids $\alpha \in \{ w, n \}$. \-Thus each component $\kappa \{ w, a \}$ can be present in each phase. \-Using the standard multiphase \-Darcy approach a mass balance equation is solved\-: \begin{eqnarray*} && \phi \frac{\partial (\sum_\alpha \varrho_\alpha X_\alpha^\kappa S_\alpha )}{\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha X_\alpha^\kappa \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} (\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g}) \right\}\\ &-& \sum_\alpha \text{div} \left\{{\bf D}_{\alpha, pm}^\kappa \varrho_{\alpha} \text{grad}\, X^\kappa_{\alpha} \right\} - \sum_\alpha q_\alpha^\kappa = 0 \qquad \kappa \in \{w, a\} \, , \alpha \in \{w, n\} \end{eqnarray*} \-For the energy balance, local thermal equilibrium is assumed which results in one energy conservation equation for the porous solid matrix and the fluids\-: \begin{eqnarray*} && \phi \frac{\partial \left( \sum_\alpha \varrho_\alpha u_\alpha S_\alpha \right)}{\partial t} + \left( 1 - \phi \right) \frac{\partial (\varrho_s c_s T)}{\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha h_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mathbf{K} \left( \text{grad}\, p_\alpha - \varrho_\alpha \mathbf{g} \right) \right\} \\ &-& \text{div} \left( \lambda_{pm} \text{grad} \, T \right) - q^h = 0 \qquad \alpha \in \{w, n\} \end{eqnarray*}
-Adaption of the BOX scheme to the non-\/isothermal two-\/phase two-\/component flow model. This model implements a non-\/isothermal two-\/phase flow of two compressible and partly miscible fluids $\alpha \in \{ w, n \}$. Thus each component $\kappa \{ w, a \}$ can be present in each phase. Using the standard multiphase Darcy approach a mass balance equation is solved: \begin{eqnarray*} && \phi \frac{\partial (\sum_\alpha \varrho_\alpha X_\alpha^\kappa S_\alpha )}{\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha X_\alpha^\kappa \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} (\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g}) \right\}\\ &-& \sum_\alpha \text{div} \left\{{\bf D_{\alpha, pm}^\kappa} \varrho_{\alpha} \text{grad}\, X^\kappa_{\alpha} \right\} - \sum_\alpha q_\alpha^\kappa = 0 \qquad \kappa \in \{w, a\} \, , \alpha \in \{w, n\} \end{eqnarray*} For the energy balance, local thermal equilibrium is assumed which results in one energy conservation equation for the porous solid matrix and the fluids: \begin{eqnarray*} && \phi \frac{\partial \left( \sum_\alpha \varrho_\alpha u_\alpha S_\alpha \right)}{\partial t} + \left( 1 - \phi \right) \frac{\partial (\varrho_s c_s T)}{\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha h_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mathbf{K} \left( \text{grad}\, p_\alpha - \varrho_\alpha \mathbf{g} \right) \right\} \\ &-& \text{div} \left( \lambda_{pm} \text{grad} \, T \right) - q^h = 0 \qquad \alpha \in \{w, n\} \end{eqnarray*}
+\-This is discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit \-Euler method as temporal discretization.
-This is discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit Euler method as temporal discretization.
-
-By using constitutive relations for the capillary pressure $p_c = p_n - p_w$ and relative permeability $k_{r\alpha}$ and taking advantage of the fact that $S_w + S_n = 1$ and $X^\kappa_w + X^\kappa_n = 1$, the number of unknowns can be reduced to two. If both phases are present the primary variables are, like in the nonisothermal two-\/phase model, either $p_w$, $S_n$ and temperature or $p_n$, $S_w$ and temperature. The formulation which ought to be used can be specified by setting the {\ttfamily Formulation} property to either {\ttfamily TwoPTwoIndices::pWsN} or {\ttfamily TwoPTwoCIndices::pNsW}. By default, the model uses $p_w$ and $S_n$. In case that only one phase (nonwetting or wetting phase) is present the second primary variable represents a mass fraction. The correct assignment of the second primary variable is performed by a phase state dependent primary variable switch. The phase state is stored for all nodes of the system. The following cases can be distinguished:
+\-By using constitutive relations for the capillary pressure $p_c = p_n - p_w$ and relative permeability $k_{r\alpha}$ and taking advantage of the fact that $S_w + S_n = 1$ and $X^\kappa_w + X^\kappa_n = 1$, the number of unknowns can be reduced to two. \-If both phases are present the primary variables are, like in the nonisothermal two-\/phase model, either $p_w$, $S_n$ and temperature or $p_n$, $S_w$ and temperature. \-The formulation which ought to be used can be specified by setting the {\ttfamily \-Formulation} property to either {\ttfamily \-Two\-P\-Two\-Indices\-::p\-Ws\-N} or {\ttfamily \-Two\-P\-Two\-C\-Indices\-::p\-Ns\-W}. \-By default, the model uses $p_w$ and $S_n$. \-In case that only one phase (nonwetting or wetting phase) is present the second primary variable represents a mass fraction. \-The correct assignment of the second primary variable is performed by a phase state dependent primary variable switch. \-The phase state is stored for all nodes of the system. \-The following cases can be distinguished\-:
\begin{itemize}
-\item Both phases are present: The saturation is used (either $S_n$ or $S_w$, dependent on the chosen formulation).
-\item Only wetting phase is present: The mass fraction of air in the wetting phase $X^a_w$ is used.
-\item Only non-\/wetting phase is present: The mass fraction of water in the non-\/wetting phase, $X^w_n$, is used.
+\item \-Both phases are present\-: \-The saturation is used (either $S_n$ or $S_w$, dependent on the chosen formulation).
+\item \-Only wetting phase is present\-: \-The mass fraction of air in the wetting phase $X^a_w$ is used.
+\item \-Only non-\/wetting phase is present\-: \-The mass fraction of water in the non-\/wetting phase, $X^w_n$, is used.
\end{itemize}
diff --git a/doc/handbook/ModelDescriptions/2pboxmodel.tex b/doc/handbook/ModelDescriptions/2pboxmodel.tex
index 9f77ea54633119ea97d769ccee6b51b508a2a9da..0bde9f6aa76bab1895e35bbdae67c1f232aac13e 100644
--- a/doc/handbook/ModelDescriptions/2pboxmodel.tex
+++ b/doc/handbook/ModelDescriptions/2pboxmodel.tex
@@ -4,10 +4,9 @@
% file instead!! %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\-This model implements two-\/phase flow of two completely immiscible fluids $\alpha \in \{ w, n \}$ using a standard multiphase \-Darcy approach as the equation for the conservation of momentum\-: \[ v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \textbf{K} \left(\text{grad}\, p_\alpha - \varrho_{\alpha} {\textbf g} \right) \]
-Adaption of the BOX scheme to the twophase flow model. This model implements two-\/phase flow of two completely immiscible fluids $\alpha \in \{ w, n \}$ using a standard multiphase Darcy approach as the equation for the conservation of momentum: \[ v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left(\text{grad} p_\alpha - \varrho_{\alpha} \mbox{\bf g} \right) \]
+\-By inserting this into the equation for the conservation of the phase mass, one gets \[ \phi \frac{\partial \varrho_\alpha S_\alpha}{\partial t} - \text{div} \left\{ \varrho_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left(\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g} \right) \right\} - q_\alpha = 0 \;, \] discretized by a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit \-Euler method as time discretization.
-By inserting this into the equation for the conservation of the phase mass, one gets \[ \phi \frac{\partial \varrho_\alpha S_\alpha}{\partial t} - \text{div} \left\{ \varrho_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left(\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g} \right) \right\} - q_\alpha = 0 \;, \] discretized by a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit Euler method as time discretization.
-
-By using constitutive relations for the capillary pressure $p_c = p_n - p_w$ and relative permeability $k_{r\alpha}$ and taking advantage of the fact that $S_w + S_n = 1$, the number of unknowns can be reduced to two. Currently the model supports choosing either $p_w$ and $S_n$ or $p_n$ and $S_w$ as primary variables. The formulation which ought to be used can be specified by setting the {\ttfamily Formulation} property to either {\ttfamily TwoPCommonIndices::pWsN} or {\ttfamily TwoPCommonIndices::pNsW}. By default, the model uses $p_w$ and $S_n$.
+\-By using constitutive relations for the capillary pressure $p_c = p_n - p_w$ and relative permeability $k_{r\alpha}$ and taking advantage of the fact that $S_w + S_n = 1$, the number of unknowns can be reduced to two. \-Currently the model supports choosing either $p_w$ and $S_n$ or $p_n$ and $S_w$ as primary variables. \-The formulation which ought to be used can be specified by setting the {\ttfamily \-Formulation} property to either {\ttfamily \-Two\-P\-Common\-Indices\-::p\-Ws\-N} or {\ttfamily \-Two\-P\-Common\-Indices\-::p\-Ns\-W}. \-By default, the model uses $p_w$ and $S_n$.
diff --git a/doc/handbook/ModelDescriptions/2pniboxmodel.tex b/doc/handbook/ModelDescriptions/2pniboxmodel.tex
index a04d33e5d862b904b2427bc35e030846c6fd79dc..1de0cb0b4ee299a08cc28799a984a6c04a2fe2af 100644
--- a/doc/handbook/ModelDescriptions/2pniboxmodel.tex
+++ b/doc/handbook/ModelDescriptions/2pniboxmodel.tex
@@ -4,10 +4,9 @@
% file instead!! %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\-This model implements a non-\/isothermal two-\/phase flow of two completely immiscible fluids $\alpha \in \{ w, n \}$. \-Using the standard multiphase \-Darcy approach, the mass conservation equations for both phases can be described as follows\-: \begin{eqnarray*} && \phi \frac{\partial (\varrho_\alpha S_\alpha )}{\partial t} - \text{div} \left\{ \varrho_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} (\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g}) \right\} - q_\alpha^\kappa = 0 \qquad \alpha \in \{w, n\} \end{eqnarray*} \-For the energy balance, local thermal equilibrium is assumed which results in one energy conservation equation for the porous solid matrix and the fluids\-: \begin{eqnarray*} && \phi \frac{\partial \left( \sum_\alpha \varrho_\alpha u_\alpha S_\alpha \right)}{\partial t} + \left( 1 - \phi \right) \frac{\partial (\varrho_s c_s T)}{\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha h_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left( \text{grad} \, p_\alpha - \varrho_\alpha \mbox{\bf g} \right) \right\} \\ &-& \text{div} \left( \lambda_{pm} \text{grad} \, T \right) - q^h = 0, \qquad \alpha \in \{w, n\}. \end{eqnarray*}
-Adaption of the BOX scheme to the non-\/isothermal twophase flow model. This model implements a non-\/isothermal two-\/phase flow of two completely immiscible fluids $\alpha \in \{ w, n \}$. Using the standard multiphase Darcy approach, the mass conservation equations for both phases can be described as follows: \begin{eqnarray*} && \phi \frac{\partial (\varrho_\alpha S_\alpha )}{\partial t} - \text{div} \left\{ \varrho_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} (\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g}) \right\} - q_\alpha^\kappa = 0 \qquad \alpha \in \{w, n\} \end{eqnarray*} For the energy balance, local thermal equilibrium is assumed which results in one energy conservation equation for the porous solid matrix and the fluids: \begin{eqnarray*} && \phi \frac{\partial \left( \sum_\alpha \varrho_\alpha u_\alpha S_\alpha \right)}{\partial t} + \left( 1 - \phi \right) \frac{\partial (\varrho_s c_s T)}{\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha h_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left( \text{grad} \, p_\alpha - \varrho_\alpha \mbox{\bf g} \right) \right\} \\ &-& \text{div} \left( \lambda_{pm} \text{grad} \, T \right) - q^h = 0, \qquad \alpha \in \{w, n\}. \end{eqnarray*}
+\-The equations are discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit \-Euler method as time discretization.
-The equations are discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit Euler method as time discretization.
-
-Currently the model supports choosing either $p_w$, $S_n$ and $T$ or $p_n$, $S_w$ and $T$ as primary variables. The formulation which ought to be used can be specified by setting the {\ttfamily Formulation} property to either {\ttfamily TwoPNIIndices::pWsN} or {\ttfamily TwoPIndices::pNsW}. By default, the model uses $p_w$, $S_n$ and $T$.
+\-Currently the model supports choosing either $p_w$, $S_n$ and $T$ or $p_n$, $S_w$ and $T$ as primary variables. \-The formulation which ought to be used can be specified by setting the {\ttfamily \-Formulation} property to either {\ttfamily \-Two\-P\-N\-I\-Indices\-::p\-Ws\-N} or {\ttfamily \-Two\-P\-Indices\-::p\-Ns\-W}. \-By default, the model uses $p_w$, $S_n$ and $T$.
diff --git a/doc/handbook/ModelDescriptions/3p3cboxmodel.tex b/doc/handbook/ModelDescriptions/3p3cboxmodel.tex
new file mode 100644
index 0000000000000000000000000000000000000000..5bd6b88d77c698c649281eeb7e95a73f08ce1550
--- /dev/null
+++ b/doc/handbook/ModelDescriptions/3p3cboxmodel.tex
@@ -0,0 +1,26 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This file has been autogenerated from the LaTeX part of the %
+% doxygen documentation; DO NOT EDIT IT! Change the model's .hh %
+% file instead!! %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+\-This model implements three-\/phase three-\/component flow of three fluid phases $\alpha \in \{ water, gas, NAPL \}$ each composed of up to three components $\kappa \in \{ water, air, contaminant \}$. \-The standard multiphase \-Darcy approach is used as the equation for the conservation of momentum\-: \[ v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left(\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g} \right) \]
+
+\-By inserting this into the equations for the conservation of the components, one transport equation for each component is obtained as \begin{eqnarray} && \phi \frac{\partial (\sum_\alpha \varrho_{\text{mol}, \alpha} x_\alpha^\kappa S_\alpha )}{\partial t} - \sum\limits_\alpha \text{div} \left\{ \frac{k_{r\alpha}}{\mu_\alpha} \varrho_{\text{mol}, \alpha} x_\alpha^\kappa \mbox{\bf K} (\text{grad}\, p_\alpha - \varrho_{\text{mass}, \alpha} \mbox{\bf g}) \right\} \nonumber \\ \nonumber \\ && - \sum\limits_\alpha \text{div} \left\{ D_{pm}^\kappa \varrho_{\text{mol}, \alpha } \text{grad}\, x_\alpha^\kappa \right\} - q^\kappa = 0 \qquad \forall \kappa , \; \forall \alpha \end{eqnarray}
+
+\-Note that these balance equations are molar.
+
+\-The equations are discretized using a fully-\/coupled vertex centered finite volume (\-B\-O\-X) scheme as spatial scheme and the implicit \-Euler method as temporal discretization.
+
+\-The model uses commonly applied auxiliary conditions like $S_w + S_n + S_g = 1$ for the saturations and $x^w_\alpha + x^a_\alpha + x^c_\alpha = 1$ for the mole fractions. \-Furthermore, the phase pressures are related to each other via capillary pressures between the fluid phases, which are functions of the saturation, e.\-g. according to the approach of \-Parker et al.
+
+\-The used primary variables are dependent on the locally present fluid phases \-An adaptive primary variable switch is included. \-The phase state is stored for all nodes of the system. \-The following cases can be distinguished\-:
+\begin{itemize}
+\item \-All three phases are present\-: \-Primary variables are two saturations ( $S_w$ and $S_n$, and a pressure, in this case $p_g$.
+\item \-Only the water phase is present\-: \-Primary variables are now the mole fractions of air and contaminant in the water phase ( $x_w^a$ and $x_w^c$), as well as the gas pressure, which is, of course, in a case where only the water phase is present, just the same as the water pressure.
+\item \-Gas and \-N\-A\-P\-L phases are present\-: \-Primary variables ( $S_n$, $x_g^w$, $p_g$).
+\item \-Water and \-N\-A\-P\-L phases are present\-: \-Primary variables ( $S_n$, $x_w^a$, $p_g$).
+\item \-Only gas phase is present\-: \-Primary variables ( $x_g^w$, $x_g^c$, $p_g$).
+\item \-Water and gas phases are present\-: \-Primary variables ( $S_w$, $x_w^g$, $p_g$).
+\end{itemize}
+
diff --git a/doc/handbook/ModelDescriptions/3p3cniboxmodel.tex b/doc/handbook/ModelDescriptions/3p3cniboxmodel.tex
new file mode 100644
index 0000000000000000000000000000000000000000..06cbf16327230451f95c06bb2229ebe75ef4fa2b
--- /dev/null
+++ b/doc/handbook/ModelDescriptions/3p3cniboxmodel.tex
@@ -0,0 +1,26 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This file has been autogenerated from the LaTeX part of the %
+% doxygen documentation; DO NOT EDIT IT! Change the model's .hh %
+% file instead!! %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+\-This model implements three-\/phase three-\/component flow of three fluid phases $\alpha \in \{ water, gas, NAPL \}$ each composed of up to three components $\kappa \in \{ water, air, contaminant \}$. \-The standard multiphase \-Darcy approach is used as the equation for the conservation of momentum\-: \[ v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \mbox{\bf K} \left(\text{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g} \right) \]
+
+\-By inserting this into the equations for the conservation of the components, one transport equation for each component is obtained as \begin{eqnarray*} && \phi \frac{\partial (\sum_\alpha \varrho_{\text{mol}, \alpha} x_\alpha^\kappa S_\alpha )}{\partial t} - \sum\limits_\alpha \text{div} \left\{ \frac{k_{r\alpha}}{\mu_\alpha} \varrho_{\text{mol}, \alpha} x_\alpha^\kappa \mbox{\bf K} (\text{grad}\; p_\alpha - \varrho_{\text{mass}, \alpha} \mbox{\bf g}) \right\} \nonumber \\ \nonumber \\ && - \sum\limits_\alpha \text{div} \left\{ D_{pm}^\kappa \varrho_{\text{mol}, \alpha } \text{grad} \; x_\alpha^\kappa \right\} - q^\kappa = 0 \qquad \forall \kappa , \; \forall \alpha \end{eqnarray*}
+
+\-Note that these balance equations above are molar. \-In addition to that, a single balance of thermal energy is formulated for the fluid-\/filled porous medium under the assumption of local thermal equilibrium \begin{eqnarray*} && \phi \frac{\partial \left( \sum_\alpha \varrho_\alpha u_\alpha S_\alpha \right)}{\partial t} + \left( 1 - \phi \right) \frac{\partial (\varrho_s c_s T)}{\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha h_\alpha \frac{k_{r\alpha}}{\mu_\alpha} \mathbf{K} \left( \text{grad}\, p_\alpha - \varrho_\alpha \mathbf{g} \right) \right\} \\ &-& \text{div} \left( \lambda_{pm} \text{grad} \, T \right) - q^h = 0 \qquad \alpha \in \{w, n, g\} \end{eqnarray*}
+
+\-The equations are discretized using a fully-\/coupled vertex centered finite volume (\-B\-O\-X) scheme as spatial scheme and the implicit \-Euler method as temporal discretization.
+
+\-The model uses commonly applied auxiliary conditions like $S_w + S_n + S_g = 1$ for the saturations and $x^w_\alpha + x^a_\alpha + x^c_\alpha = 1$ for the mole fractions. \-Furthermore, the phase pressures are related to each other via capillary pressures between the fluid phases, which are functions of the saturation, e.\-g. according to the approach of \-Parker et al.
+
+\-The used primary variables are dependent on the locally present fluid phases \-An adaptive primary variable switch is included. \-The phase state is stored for all nodes of the system. \-The following cases can be distinguished\-:
+\begin{itemize}
+\item \-All three phases are present\-: \-Primary variables are two saturations ( $S_w$ and $S_n$, a pressure, in this case $p_g$, and the temperature $T$.
+\item \-Only the water phase is present\-: \-Primary variables are now the mole fractions of air and contaminant in the water phase ( $x_w^a$ and $x_w^c$), as well as temperature and the gas pressure, which is, of course, in a case where only the water phase is present, just the same as the water pressure.
+\item \-Gas and \-N\-A\-P\-L phases are present\-: \-Primary variables ( $S_n$, $x_g^w$, $p_g$, $T$).
+\item \-Water and \-N\-A\-P\-L phases are present\-: \-Primary variables ( $S_n$, $x_w^a$, $p_g$, $T$).
+\item \-Only gas phase is present\-: \-Primary variables ( $x_g^w$, $x_g^c$, $p_g$, $T$).
+\item \-Water and gas phases are present\-: \-Primary variables ( $S_w$, $x_w^g$, $p_g$, $T$).
+\end{itemize}
+
diff --git a/doc/handbook/ModelDescriptions/richardsboxmodel.tex b/doc/handbook/ModelDescriptions/richardsboxmodel.tex
index bdf9c47e43128e7679413459f79c6858b9a25957..18195eff89b98be8e99d8598f0f4aa27d60488e7 100644
--- a/doc/handbook/ModelDescriptions/richardsboxmodel.tex
+++ b/doc/handbook/ModelDescriptions/richardsboxmodel.tex
@@ -4,10 +4,9 @@
% file instead!! %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\-In the unsaturated zone, \-Richards' equation is frequently used to calculate the water distribution above the groundwater level. \-It can be derived from the twophase equations, i.\-e. \[ \frac{\partial\;\phi S_\alpha \rho_\alpha}{\partial t} - \text{div} \left\{ \rho_\alpha \frac{k_{r\alpha}}{\mu_\alpha}\; {\textbf K} \text{grad}\left[ p_\alpha - g\rho_\alpha \right] \right\} = q_\alpha, \] where $\alpha \in \{w, n\}$ is the fluid phase, $\rho_\alpha$ is the fluid density, $S_\alpha$ is the fluid saturation, $\phi$ is the porosity of the soil, $k_{r\alpha}$ is the relative permeability for the fluid, $\mu_\alpha$ is the fluid's dynamic viscosity, $K$ is the intrinsic permeability, $p_\alpha$ is the fluid pressure and $g$ is the potential of the gravity field.
-This model implements a variant of the Richards equation for quasi-\/twophase flow. In the unsaturated zone, Richards' equation is frequently used to calculate the water distribution above the groundwater level. It can be derived from the twophase equations, i.e. \[ \frac{\partial\;\phi S_\alpha \rho_\alpha}{\partial t} - \mathbf{div} \left\{ \rho_\alpha \frac{k_{r\alpha}}{\mu_\alpha}\;K \mathbf{grad}\left[ p_\alpha - g\rho_\alpha \right] \right\} = q_\alpha, \] where $\alpha \in \{w, n\}$ is the fluid phase, $\rho_\alpha$ is the fluid density, $S_\alpha$ is the fluid saturation, $\phi$ is the porosity of the soil, $k_{r\alpha}$ is the relative permeability for the fluid, $\mu_\alpha$ is the fluid's dynamic viscosity, $K$ is the intrinsic permeability, $p_\alpha$ is the fluid pressure and $g$ is the potential of the gravity field.
+\-In contrast to the full twophase model, the \-Richards model assumes gas as the non-\/wetting fluid and that it exhibits a much lower viscosity than the (liquid) wetting phase. (\-For example at atmospheric pressure and at room temperature, the viscosity of air is only about $1\%$ of the viscosity of liquid water.) \-As a consequence, the $\frac{k_{r\alpha}}{\mu_\alpha}$ term typically is much larger for the gas phase than for the wetting phase. \-For this reason, the \-Richards model assumes that $\frac{k_{rn}}{\mu_n}$ tends to infinity. \-This implies that the pressure of the gas phase is equivalent to a static pressure and can thus be specified externally and that therefore, mass conservation only needs to be considered for the wetting phase.
-In contrast to the full twophase model, the Richards model assumes gas as the non-\/wetting fluid and that it exhibits a much lower viscosity than the (liquid) wetting phase. (For example at atmospheric pressure and at room temperature, the viscosity of air is only about $1\%$ of the viscosity of liquid water.) As a consequence, the $\frac{k_{r\alpha}}{\mu_\alpha}$ term typically is much larger for the gas phase than for the wetting phase. For this reason, the Richards model assumes that $\frac{k_{rn}}{\mu_n}$ tends to infinity. This implies that the pressure of the gas phase is equivalent to a static pressure and can thus be specified externally and that therefore, mass conservation only needs to be considered for the wetting phase.
-
-The model thus choses the absolute pressure of the wetting phase $p_w$ as its only primary variable. The wetting phase saturation is calculated using the inverse of the capillary pressure, i.e. \[ S_w = p_c^{-1}(p_n - p_w) \] holds, where $p_n$ is a given reference pressure. Nota bene that the last step is assumes that the capillary pressure-\/saturation curve can be inverted uniquely, so it is not possible to set the capillary pressure to zero when using the Richards model!
+\-The model thus choses the absolute pressure of the wetting phase $p_w$ as its only primary variable. \-The wetting phase saturation is calculated using the inverse of the capillary pressure, i.\-e. \[ S_w = p_c^{-1}(p_n - p_w) \] holds, where $p_n$ is a given reference pressure. \-Nota bene that the last step is assumes that the capillary pressure-\/saturation curve can be inverted uniquely, so it is not possible to set the capillary pressure to zero when using the \-Richards model!
diff --git a/doc/handbook/ModelDescriptions/stokes2cmodel.tex b/doc/handbook/ModelDescriptions/stokes2cmodel.tex
new file mode 100644
index 0000000000000000000000000000000000000000..9e9dd9d58626edfd06ef18508f60ff52754e7bbd
--- /dev/null
+++ b/doc/handbook/ModelDescriptions/stokes2cmodel.tex
@@ -0,0 +1,16 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This file has been autogenerated from the LaTeX part of the %
+% doxygen documentation; DO NOT EDIT IT! Change the model's .hh %
+% file instead!! %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+\-This model implements an isothermal two-\/component \-Stokes flow of a fluid solving a momentum balance, a mass balance and a conservation equation for one component.
+
+\-Momentum \-Balance\-: \[ \frac{\partial \left(\varrho_g {\boldsymbol{v}}_g\right)}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} \left(p_g {\bf {I}} - \mu_g \left(\boldsymbol{\nabla} \boldsymbol{v}_g + \boldsymbol{\nabla} \boldsymbol{v}_g^T\right)\right) - \varrho_g {\bf g} = 0, \]
+
+\-Mass balance equation\-: \[ \frac{\partial \varrho_g}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\varrho_g {\boldsymbol{v}}_g\right) - q_g = 0 \]
+
+\hyperlink{a00047}{\-Component} mass balance equation\-: \[ \frac{\partial \left(\varrho_g X_g^\kappa\right)}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} \left( \varrho_g {\boldsymbol{v}}_g X_g^\kappa - D^\kappa_g \varrho_g \boldsymbol{\nabla} X_g^\kappa \right) - q_g^\kappa = 0 \]
+
+\-This is discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit \-Euler method as temporal discretization.
+
diff --git a/doc/handbook/ModelDescriptions/stokes2cnimodel.tex b/doc/handbook/ModelDescriptions/stokes2cnimodel.tex
new file mode 100644
index 0000000000000000000000000000000000000000..ef8c8cc4ad5b1e0de69c99293b7a916b92e67925
--- /dev/null
+++ b/doc/handbook/ModelDescriptions/stokes2cnimodel.tex
@@ -0,0 +1,18 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This file has been autogenerated from the LaTeX part of the %
+% doxygen documentation; DO NOT EDIT IT! Change the model's .hh %
+% file instead!! %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+\-This model implements a non-\/isothermal two-\/component \-Stokes flow of a fluid solving a momentum balance, a mass balance, a conservation equation for one component, and one balance equation for the energy.
+
+\-Momentum \-Balance\-: \[ \frac{\partial \left(\varrho_g {\boldsymbol{v}}_g\right)}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} \left(p_g {\bf {I}} - \mu_g \left(\boldsymbol{\nabla} \boldsymbol{v}_g + \boldsymbol{\nabla} \boldsymbol{v}_g^T\right)\right) - \varrho_g {\bf g} = 0, \]
+
+\-Mass balance equation\-: \[ \frac{\partial \varrho_g}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\varrho_g {\boldsymbol{v}}_g\right) - q_g = 0 \]
+
+\hyperlink{a00047}{\-Component} mass balance equation\-: \[ \frac{\partial \left(\varrho_g X_g^\kappa\right)}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} \left( \varrho_g {\boldsymbol{v}}_g X_g^\kappa - D^\kappa_g \varrho_g \boldsymbol{\nabla} X_g^\kappa \right) - q_g^\kappa = 0 \]
+
+\-Energy balance equation\-: \[ \frac{\partial (\varrho_g u_g)}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} \varrho_g h_g {\boldsymbol{v}}_g - \lambda_g \boldsymbol{\nabla} T - q_T = 0 \]
+
+\-This is discretized using a fully-\/coupled vertex centered finite volume (box) scheme as spatial and the implicit \-Euler method as temporal discretization.
+
diff --git a/doc/handbook/ModelDescriptions/stokesmodel.tex b/doc/handbook/ModelDescriptions/stokesmodel.tex
new file mode 100644
index 0000000000000000000000000000000000000000..b855afa3b5eac9171219e48ddfb816714b1d9b2e
--- /dev/null
+++ b/doc/handbook/ModelDescriptions/stokesmodel.tex
@@ -0,0 +1,12 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This file has been autogenerated from the LaTeX part of the %
+% doxygen documentation; DO NOT EDIT IT! Change the model's .hh %
+% file instead!! %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+\-This model implements laminar \-Stokes flow of a single fluid, solving a momentum balance\-: \[ \frac{\partial \left(\varrho_g {\boldsymbol{v}}_g\right)}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} \left(p_g {\bf {I}} - \mu_g \left(\boldsymbol{\nabla} \boldsymbol{v}_g + \boldsymbol{\nabla} \boldsymbol{v}_g^T\right)\right) - \varrho_g {\bf g} = 0, \]
+
+and the mass balance equation\-: \[ \frac{\partial \varrho_g}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\varrho_g {\boldsymbol{v}}_g\right) - q_g = 0 \]
+
+\-This is discretized by a fully-\/coupled vertex-\/centered finite volume (box) scheme in space and by the implicit \-Euler method in time.
+
diff --git a/doc/handbook/models.tex b/doc/handbook/models.tex
index c117929653b1f18abcb72a96faa139fff1565167..a16f19f4b932fda7a29801f7fc3f785286a6b3f0 100644
--- a/doc/handbook/models.tex
+++ b/doc/handbook/models.tex
@@ -147,7 +147,8 @@ fractions, temperature, pressures, etc.
\section{Available models}
The following description of the available models is automatically extracted
-from the Doxygen documentation. \textbf{TODO}: Unify notation.
+from the Doxygen documentation.
+% \textbf{TODO}: Unify notation.
\subsection{Fully coupled models}
@@ -157,7 +158,7 @@ from the Doxygen documentation. \textbf{TODO}: Unify notation.
\subsubsection{The single-phase, two-component model: OnePTwoCBoxModel}
\input{ModelDescriptions/1p2cboxmodel}
-\subsubsection{The two-phase model using Richards' assumption: RichardsBoxModel}
+\subsubsection{The two-phase model using the \textsc{Richards} assumption: RichardsBoxModel}
\input{ModelDescriptions/richardsboxmodel}
\subsubsection{The two-phase model: TwoPBoxModel}
@@ -172,6 +173,20 @@ from the Doxygen documentation. \textbf{TODO}: Unify notation.
\subsubsection{The non-isothermal two-phase, two-component model: TwoPTwoCNIBoxModel}
\input{ModelDescriptions/2p2cniboxmodel}
+\subsubsection{The three-phase, three-component model: ThreePThreeCNIBoxModel}
+\input{ModelDescriptions/3p3cboxmodel}
+
+\subsubsection{The non-isothermal three-phase, three-component model: ThreePThreeCNIBoxModel}
+\input{ModelDescriptions/3p3cniboxmodel}
+
+\subsubsection{The \textsc{Stokes} model: StokesModel}
+\input{ModelDescriptions/stokesmodel}
+
+\subsubsection{The isothermal two-component \textsc{Stokes} model: Stokes2cModel}
+\input{ModelDescriptions/stokes2cmodel}
+
+\subsubsection{The non-isothermal two-component \textsc{Stokes} model: Stokes2cniModel}
+\input{ModelDescriptions/stokes2cnimodel}
\subsection{Decoupled models}
@@ -181,8 +196,6 @@ from the Doxygen documentation. \textbf{TODO}: Unify notation.
\subsubsection{Decoupled Compositional Models}
\input{ModelDescriptions/decoupled2p2c}
-
-
%%% Local Variables:
%%% mode: latex
%%% TeX-master: "dumux-handbook"
diff --git a/dumux/boxmodels/2p/2pmodel.hh b/dumux/boxmodels/2p/2pmodel.hh
index c85084449da779e9f2bd868ea9e29cd5ce1d401e..82f88311a81b679dbea22d357a1e3e6785591d16 100644
--- a/dumux/boxmodels/2p/2pmodel.hh
+++ b/dumux/boxmodels/2p/2pmodel.hh
@@ -45,7 +45,7 @@ namespace Dumux
* \f$\alpha \in \{ w, n \}\f$ using a standard multiphase Darcy
* approach as the equation for the conservation of momentum:
\f[
- v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \textbf K}
+ v_\alpha = - \frac{k_{r\alpha}}{\mu_\alpha} \textbf{K}
\left(\text{grad}\, p_\alpha - \varrho_{\alpha} {\textbf g} \right)
\f]
*