3.2-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
Classes | Namespaces
porousmediumflow/2p2c/model.hh File Reference

Adaption of the fully implicit scheme to the two-phase two-component fully implicit model. More...

#include <array>
#include <dumux/common/properties.hh>
#include <dumux/porousmediumflow/2pnc/model.hh>
#include <dumux/porousmediumflow/2p/formulation.hh>
#include <dumux/porousmediumflow/nonisothermal/model.hh>
#include <dumux/porousmediumflow/nonisothermal/iofields.hh>
#include <dumux/porousmediumflow/nonequilibrium/model.hh>
#include <dumux/porousmediumflow/nonequilibrium/volumevariables.hh>
#include <dumux/material/fluidmatrixinteractions/2p/thermalconductivitysomerton.hh>
#include <dumux/material/fluidmatrixinteractions/2p/thermalconductivitysimplefluidlumping.hh>
#include "volumevariables.hh"

Go to the source code of this file.

Description

Adaption of the fully implicit scheme to the two-phase two-component fully implicit 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} \mathbf{K} \left(\textbf{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 \frac{M^\kappa}{M_\alpha} x_\alpha^\kappa S_\alpha )} {\partial t} - \sum_\alpha \text{div} \left\{ \varrho_\alpha \frac{M^\kappa}{M_\alpha} x_\alpha^\kappa \frac{k_{r\alpha}}{\mu_\alpha} \mathbf{K} (\textbf{grad}\, p_\alpha - \varrho_{\alpha} \mbox{\bf g}) \right\} \nonumber \\ \nonumber \\ &-& \sum_\alpha \text{div} \left\{ D_{\alpha,\text{pm}}^\kappa \varrho_{\alpha} \textbf{grad} X^\kappa_{\alpha} \right\} - \sum_\alpha q_\alpha^\kappa = 0 \qquad \kappa \in \{w, a\} \, , \alpha \in \{w, g\} \end{eqnarray*}

All equations are discretized using a vertex-centered finite volume (box) or cell-centered finite volume 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\) 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 Formulation property to either TwoPTwoCFormulation::pwsn or TwoPTwoCFormulation::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. The model is able to use either mole or mass fractions. The property useMoles can be set to either true or false in the problem file. Make sure that the according units are used in the problem setup. useMoles is set to true by default. Following cases can be distinguished:

Classes

struct  Dumux::Properties::TTag::TwoPTwoC
 
struct  Dumux::Properties::TTag::TwoPTwoCNI
 
struct  Dumux::Properties::BaseModelTraits< TypeTag, TTag::TwoPTwoC >
 Set the model traits property. More...
 
struct  Dumux::Properties::ModelTraits< TypeTag, TTag::TwoPTwoC >
 
struct  Dumux::Properties::VolumeVariables< TypeTag, TTag::TwoPTwoC >
 Use the 2p2c VolumeVariables. More...
 
struct  Dumux::Properties::UseConstraintSolver< TypeTag, TTag::TwoPTwoC >
 Determines whether the constraint solver is used. More...
 
struct  Dumux::Properties::ModelTraits< TypeTag, TTag::TwoPTwoCNI >
 Set the non-isothermal model traits. More...
 
struct  Dumux::Properties::VolumeVariables< TypeTag, TTag::TwoPTwoCNI >
 Set the volume variables property. More...
 
struct  Dumux::Properties::IOFields< TypeTag, TTag::TwoPTwoCNI >
 Set non-isothermal output fields. More...
 
struct  Dumux::Properties::ThermalConductivityModel< TypeTag, TTag::TwoPTwoCNI >
 Somerton is used as default model to compute the effective thermal heat conductivity. More...
 
struct  Dumux::TwoPTwoCUnconstrainedModelTraits< TwoPTwoCModelTraits >
 
struct  Dumux::Properties::TTag::TwoPTwoCNonEquil
 
struct  Dumux::Properties::EquilibriumLocalResidual< TypeTag, TTag::TwoPTwoCNonEquil >
 
struct  Dumux::Properties::EquilibriumIOFields< TypeTag, TTag::TwoPTwoCNonEquil >
 Set the vtk output fields specific to this model. More...
 
struct  Dumux::Properties::ModelTraits< TypeTag, TTag::TwoPTwoCNonEquil >
 
struct  Dumux::Properties::EquilibriumModelTraits< TypeTag, TTag::TwoPTwoCNonEquil >
 Set equilibrium model traits. More...
 
struct  Dumux::Properties::ThermalConductivityModel< TypeTag, TTag::TwoPTwoCNonEquil >
 
struct  Dumux::Properties::VolumeVariables< TypeTag, TTag::TwoPTwoCNonEquil >
 Use the nonequilibrium volume variables together with the 2p2c vol vars. More...
 
struct  Dumux::Properties::TTag::TwoPTwoCNINonEquil
 
struct  Dumux::Properties::ModelTraits< TypeTag, TTag::TwoPTwoCNINonEquil >
 Set the non-isothermal model traits with the nonequilibrium model traits as isothermal traits. More...
 
struct  Dumux::Properties::EquilibriumIOFields< TypeTag, TTag::TwoPTwoCNINonEquil >
 Set the equilibrium IO fields which are in that case the nonisothermal io fields. More...
 
struct  Dumux::Properties::VolumeVariables< TypeTag, TTag::TwoPTwoCNINonEquil >
 Use the nonequilibrium volume variables together with the 2p2c vol vars. More...
 
struct  Dumux::Properties::ThermalConductivityModel< TypeTag, TTag::TwoPTwoCNINonEquil >
 Somerton is used as default model to compute the effective thermal heat conductivity. More...
 

Namespaces

namespace  Dumux
 
namespace  Dumux::Properties
 
namespace  Dumux::Properties::TTag
 Type tag for numeric models.
 
Include dependency graph for porousmediumflow/2p2c/model.hh: