3.4
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
Public Types | Static Public Member Functions | Static Public Attributes | List of all members
Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false > Class Template Reference

Specialization of the CCTpfaForchheimersLaw grids where dim=dimWorld. More...

#include <dumux/flux/cctpfa/forchheimerslaw.hh>

Description

template<class ScalarType, class GridGeometry, class UpwindScheme>
class Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >

Specialization of the CCTpfaForchheimersLaw grids where dim=dimWorld.

Public Types

using Scalar = ScalarType
 state the scalar type of the law More...
 
using Cache = TpfaForchheimersLawCache< ThisType, GridGeometry >
 state the type for the corresponding cache More...
 

Static Public Member Functions

template<class Problem , class ElementVolumeVariables , class ElementFluxVarsCache >
static Scalar flux (const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, int phaseIdx, const ElementFluxVarsCache &elemFluxVarsCache)
 Compute the advective flux of a phase across the given sub-control volume face uing the Forchheimer equation. More...
 
template<class Problem , class ElementVolumeVariables >
static Scalar calculateTransmissibility (const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf)
 
template<class Problem , class ElementVolumeVariables , bool scalarPerm = std::is_same<typename Problem::SpatialParams::PermeabilityType, Scalar>::value, std::enable_if_t< scalarPerm, int > = 0>
static DimWorldMatrix calculateHarmonicMeanSqrtPermeability (const Problem &problem, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf)
 Returns the harmonic mean of \(\sqrt{K_0}\) and \(\sqrt{K_1}\). More...
 
template<class Problem , class ElementVolumeVariables , bool scalarPerm = std::is_same<typename Problem::SpatialParams::PermeabilityType, Scalar>::value, std::enable_if_t<!scalarPerm, int > = 0>
static DimWorldMatrix calculateHarmonicMeanSqrtPermeability (const Problem &problem, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf)
 Returns the harmonic mean of \(\sqrt{\mathbf{K_0}}\) and \(\sqrt{\mathbf{K_1}}\). More...
 

Static Public Attributes

static const DiscretizationMethod discMethod = DiscretizationMethod::cctpfa
 state the discretization method this implementation belongs to More...
 

Member Typedef Documentation

◆ Cache

template<class ScalarType , class GridGeometry , class UpwindScheme >
using Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >::Cache = TpfaForchheimersLawCache<ThisType, GridGeometry>

state the type for the corresponding cache

◆ Scalar

template<class ScalarType , class GridGeometry , class UpwindScheme >
using Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >::Scalar = ScalarType

state the scalar type of the law

Member Function Documentation

◆ calculateHarmonicMeanSqrtPermeability() [1/2]

template<class ScalarType , class GridGeometry , class UpwindScheme >
template<class Problem , class ElementVolumeVariables , bool scalarPerm = std::is_same<typename Problem::SpatialParams::PermeabilityType, Scalar>::value, std::enable_if_t< scalarPerm, int > = 0>
static DimWorldMatrix Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >::calculateHarmonicMeanSqrtPermeability ( const Problem &  problem,
const ElementVolumeVariables &  elemVolVars,
const SubControlVolumeFace &  scvf 
)
inlinestatic

Returns the harmonic mean of \(\sqrt{K_0}\) and \(\sqrt{K_1}\).

This is a specialization for scalar-valued permeabilities which returns a tensor with identical diagonal entries.

◆ calculateHarmonicMeanSqrtPermeability() [2/2]

template<class ScalarType , class GridGeometry , class UpwindScheme >
template<class Problem , class ElementVolumeVariables , bool scalarPerm = std::is_same<typename Problem::SpatialParams::PermeabilityType, Scalar>::value, std::enable_if_t<!scalarPerm, int > = 0>
static DimWorldMatrix Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >::calculateHarmonicMeanSqrtPermeability ( const Problem &  problem,
const ElementVolumeVariables &  elemVolVars,
const SubControlVolumeFace &  scvf 
)
inlinestatic

Returns the harmonic mean of \(\sqrt{\mathbf{K_0}}\) and \(\sqrt{\mathbf{K_1}}\).

This is a specialization for tensor-valued permeabilities.

◆ calculateTransmissibility()

template<class ScalarType , class GridGeometry , class UpwindScheme >
template<class Problem , class ElementVolumeVariables >
static Scalar Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >::calculateTransmissibility ( const Problem &  problem,
const Element &  element,
const FVElementGeometry &  fvGeometry,
const ElementVolumeVariables &  elemVolVars,
const SubControlVolumeFace &  scvf 
)
inlinestatic

The flux variables cache has to be bound to an element prior to flux calculations During the binding, the transmissibility will be computed and stored using the method below.

◆ flux()

template<class ScalarType , class GridGeometry , class UpwindScheme >
template<class Problem , class ElementVolumeVariables , class ElementFluxVarsCache >
static Scalar Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >::flux ( const Problem &  problem,
const Element &  element,
const FVElementGeometry &  fvGeometry,
const ElementVolumeVariables &  elemVolVars,
const SubControlVolumeFace &  scvf,
int  phaseIdx,
const ElementFluxVarsCache &  elemFluxVarsCache 
)
inlinestatic

Compute the advective flux of a phase across the given sub-control volume face uing the Forchheimer equation.

The flux is given in N*m, and can be converted into a volume flux (m^3/s) or mass flux (kg/s) by applying an upwind scheme for the mobility or the product of density and mobility, respectively.

see e.g. Nield & Bejan: Convection in Porous Media [43]

The relative passability \( \eta_r\) is the "Forchheimer-equivalent" to the relative permeability \( k_r\). We use the same function as for \( k_r \) (VG, BC, linear) other authors use a simple power law e.g.: \(\eta_{rw} = S_w^3\)

Some rearrangements have been made to the original Forchheimer relation:

\[ \mathbf{v_\alpha} + c_F \sqrt{\mathbf{K}} \frac{\rho_\alpha}{\mu_\alpha } |\mathbf{v_\alpha}| \mathbf{v_\alpha} + \frac{k_{r \alpha}}{\mu_\alpha} \mathbf{K} \nabla \left(p_\alpha + \rho_\alpha g z \right)= 0 \]

This already includes the assumption \( k_r(S_w) = \eta_r(S_w)\):

  • \(\eta_{rw} = S_w^x\) looks very similar to e.g. Van Genuchten relative permeabilities
  • Fichot, et al. (2006), Nuclear Engineering and Design, state that several authors claim that \( k_r(S_w), \eta_r(S_w)\) can be chosen equal
  • It leads to the equation not degenerating for the case of \(S_w=1\), because I do not need to multiply with two different functions, and therefore there are terms not being zero.
  • If this assumption is not to be made: Some regularization needs to be introduced ensuring that not all terms become zero for \(S_w=1\).

This non-linear equations is solved for \(\mathbf{v_\alpha}\) using Newton's method and an analytical derivative w.r.t. \(\mathbf{v_\alpha}\).

The gradient of the Forchheimer relations looks as follows (mind that \(\sqrt{\mathbf{K}}\) is a tensor):

\[ f\left(\mathbf{v_\alpha}\right) = \left( \begin{array}{ccc} 1 & 0 &0 \\ 0 & 1 &0 \\ 0 & 0 &1 \\ \end{array} \right) + c_F \frac{\rho_\alpha}{\mu_\alpha} |\mathbf{v}_\alpha| \sqrt{\mathbf{K}} + c_F \frac{\rho_\alpha}{\mu_\alpha}\frac{1}{|\mathbf{v}_\alpha|} \sqrt{\mathbf{K}} \left( \begin{array}{ccc} v_x^2 & v_xv_y & v_xv_z \\ v_yv_x & v_{y}^2 & v_yv_z \\ v_zv_x & v_zv_y &v_{z}^2 \\ \end{array} \right) \]

Note
We restrict the use of Forchheimer's law to diagonal permeability tensors so far. This might be changed to general tensors using eigenvalue decomposition to get \(\sqrt{\mathbf{K}}\)

Member Data Documentation

◆ discMethod

template<class ScalarType , class GridGeometry , class UpwindScheme >
const DiscretizationMethod Dumux::CCTpfaForchheimersLaw< ScalarType, GridGeometry, UpwindScheme, false >::discMethod = DiscretizationMethod::cctpfa
static

state the discretization method this implementation belongs to


The documentation for this class was generated from the following file: