3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
Public Types | Static Public Member Functions | List of all members
Dumux::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables > Class Template Reference

Forchheimer's law This file contains the calculation of the Forchheimer velocity for a given Darcy velocity. More...

#include <dumux/flux/forchheimervelocity.hh>

Description

template<class Scalar, class GridGeometry, class FluxVariables>
class Dumux::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables >

Forchheimer's law This file contains the calculation of the Forchheimer velocity for a given Darcy velocity.

Public Types

using UpwindScheme = typename FluxVariables::UpwindScheme
 
using DimWorldVector = Dune::FieldVector< Scalar, dimWorld >
 
using DimWorldMatrix = Dune::FieldMatrix< Scalar, dimWorld, dimWorld >
 

Static Public Member Functions

template<class ElementVolumeVariables >
static DimWorldVector velocity (const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const int phaseIdx, const DimWorldMatrix sqrtK, const DimWorldVector darcyVelocity)
 Compute the Forchheimer velocity of a phase at the given sub-control volume face using the Forchheimer equation. 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{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...
 

Member Typedef Documentation

◆ DimWorldMatrix

template<class Scalar , class GridGeometry , class FluxVariables >
using Dumux::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables >::DimWorldMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>

◆ DimWorldVector

template<class Scalar , class GridGeometry , class FluxVariables >
using Dumux::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables >::DimWorldVector = Dune::FieldVector<Scalar, dimWorld>

◆ UpwindScheme

template<class Scalar , class GridGeometry , class FluxVariables >
using Dumux::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables >::UpwindScheme = typename FluxVariables::UpwindScheme

Member Function Documentation

◆ calculateHarmonicMeanSqrtPermeability() [1/2]

template<class Scalar , class GridGeometry , class FluxVariables >
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::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables >::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 Scalar , class GridGeometry , class FluxVariables >
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::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables >::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.

◆ velocity()

template<class Scalar , class GridGeometry , class FluxVariables >
template<class ElementVolumeVariables >
static DimWorldVector Dumux::ForchheimerVelocity< Scalar, GridGeometry, FluxVariables >::velocity ( const FVElementGeometry &  fvGeometry,
const ElementVolumeVariables &  elemVolVars,
const SubControlVolumeFace &  scvf,
const int  phaseIdx,
const DimWorldMatrix  sqrtK,
const DimWorldVector  darcyVelocity 
)
inlinestatic

Compute the Forchheimer velocity of a phase at the given sub-control volume face using the Forchheimer equation.

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

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}}\)

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