The finite volume model for the solution of the compositional pressure equation. More...
#include <dumux/porousmediumflow/2p2c/sequential/fvpressuremultiphysics.hh>
The finite volume model for the solution of the compositional pressure equation.
Provides a Finite Volume implementation for the pressure equation of a gas-liquid system with two components. An IMPES-like method is used for the sequential solution of the problem. Diffusion is neglected, capillarity can be regarded. Isothermal conditions and local thermodynamic equilibrium are assumed. Gravity is included.
\[ c_{total}\frac{\partial p}{\partial t} + \sum_{\kappa} \frac{\partial v_{total}}{\partial C^{\kappa}} \nabla \cdot \left( \sum_{\alpha} X^{\kappa}_{\alpha} \varrho_{alpha} \bf{v}_{\alpha}\right) = \sum_{\kappa} \frac{\partial v_{total}}{\partial C^{\kappa}} q^{\kappa}, \]
where \(\bf{v}_{\alpha} = - \lambda_{\alpha} \bf{K} \left(\nabla p_{\alpha} + \rho_{\alpha} \bf{g} \right) \). \( c_{total} \) represents the total compressibility, for constant porosity this yields \( - \frac{\partial V_{total}}{\partial p_{\alpha}} \), \(p_{\alpha} \) denotes the phase pressure, \( \bf{K} \) the absolute permeability, \( \lambda_{\alpha} \) the phase mobility, \( \rho_{\alpha} \) the phase density and \( \bf{g} \) the gravity constant and \( C^{\kappa} \) the total Component concentration. See paper SPE 99619 or "Analysis of a Compositional Model for Fluid Flow in Porous Media" by Chen, Qin and Ewing for derivation.
The partial derivatives of the actual fluid volume \( v_{total} \) are gained by using a secant method.
The model domain is automatically divided in a single-phase and a two-phase domain. The full 2p2c model is only evaluated within the two-phase subdomain, whereas a single-phase transport model is computed in the rest of the domain.
TypeTag | The Type Tag |
Public Member Functions | |
void | assemble (bool first) |
function which assembles the system of equations to be solved More... | |
void | get1pSource (EntryType &sourceEntry, const Element &elementI, const CellData &cellDataI) |
Assembles the source term. More... | |
void | get1pStorage (EntryType &storageEntry, const Element &elementI, CellData &cellDataI) |
Assembles the storage term for a 1p cell in a multiphysics framework. More... | |
void | get1pFlux (EntryType &entries, const Intersection &intersection, const CellData &cellDataI) |
The compositional single-phase flux in the multiphysics framework. More... | |
void | get1pFluxOnBoundary (EntryType &entries, const Intersection &intersection, const CellData &cellDataI) |
The compositional single-phase flux in the multiphysics framework. More... | |
void | initialize (bool solveTwice=false) |
void | serializeEntity (std::ostream &outstream, const Element &element) |
Function for serialization of the pressure field. More... | |
void | deserializeEntity (std::istream &instream, const Element &element) |
Function for deserialization of the pressure field. More... | |
void | updateMaterialLaws (bool postTimeStep=false) |
constitutive functions are updated once if new concentrations are calculated and stored in the variables container More... | |
void | update1pMaterialLawsInElement (const Element &elementI, CellData &cellData, bool postTimeStep) |
updates secondary variables of one single phase cell More... | |
template<class MultiWriter > | |
void | addOutputVtkFields (MultiWriter &writer) |
Write data files. More... | |
FVPressure2P2CMultiPhysics (Problem &problem) | |
Constructs a FVPressure2P2CPC object. More... | |
void | getSource (EntryType &sourceEntry, const Element &elementI, const CellData &cellDataI, const bool first) |
Assembles the source term. More... | |
void | getStorage (EntryType &storageEntry, const Element &elementI, const CellData &cellDataI, const bool first) |
Assembles the storage term. More... | |
void | getFlux (EntryType &entries, const Intersection &intersection, const CellData &cellDataI, const bool first) |
Get flux at an interface between two cells. More... | |
void | getFluxOnBoundary (EntryType &entries, const Intersection &intersection, const CellData &cellDataI, const bool first) |
Get flux on Boundary. More... | |
void | updateMaterialLawsInElement (const Element &elementI, bool postTimeStep) |
Updates secondary variables of one cell. More... | |
void | initialMaterialLaws (bool compositional) |
initializes the fluid distribution and hereby the variables container More... | |
void | initialize () |
Initialize pressure model. More... | |
void | update () |
Compositional pressure solution routine: update estimate for secants, assemble, solve. More... | |
void | volumeDerivatives (const GlobalPosition &, const Element &ep) |
Partial derivatives of the volumes w.r.t. changes in total concentration and pressure. More... | |
const Scalar | pressure (const int eIdxGlobal) const |
Public access function for the primary pressure variable. More... | |
const Matrix & | globalMatrix () |
Returns the global matrix of the last pressure solution step. More... | |
const RHSVector & | rightHandSide () |
Returns the right hand side vector of the last pressure solution step. More... | |
void | calculateVelocity () |
void | updateVelocity () |
void | setFixPressureAtIndex (Scalar pressure, int eIdxGlobal) |
Set a pressure to be fixed at a certain cell. More... | |
void | unsetFixPressureAtIndex (int eIdxGlobal) |
Reset the fixed pressure state. More... | |
void | resetFixPressureAtIndex () |
Protected Types | |
enum | { rhs = 1 , matrix = 0 } |
Indices of matrix and rhs entries. More... | |
using | ElementMapper = typename SolutionTypes::ElementMapper |
using | DataHandle = VectorExchange< ElementMapper, Dune::BlockVector< Dune::FieldVector< int, 1 > > > |
enum | { rhs = 1 , matrix = 0 } |
Indices of matrix and rhs entries. More... | |
enum | { pressEqIdx = Indices::pressureEqIdx } |
Protected Member Functions | |
void | initializeMatrix () |
Initialize the global matrix of the system of equations to solve. More... | |
void | initializeMatrixRowSize () |
Initialize the global matrix of the system of equations to solve. More... | |
void | initializeMatrixIndices () |
Initialize the global matrix of the system of equations to solve. More... | |
void | solve () |
Solves the global system of equations to get the spatial distribution of the pressure. More... | |
PressureSolution & | pressure () |
Returns the vector containing the pressure solution. More... | |
const PressureSolution & | pressure () const |
Returns the vector containing the pressure solution. More... | |
void | initializePressure () |
Initialization of the pressure solution vector: Initialization with meaningful values may. More... | |
Protected Attributes | |
Dune::BlockVector< Dune::FieldVector< int, 1 > > | nextSubdomain |
vector holding next subdomain More... | |
const GlobalPosition & | gravity_ |
Dune::Timer | timer_ |
A timer for the time spent on the multiphysics framework. More... | |
Problem & | problem_ |
bool | enableVolumeIntegral |
Enables the volume integral of the pressure equation. More... | |
bool | regulateBoundaryPermeability |
Enables regulation of permeability in the direction of a Dirichlet Boundary Condition. More... | |
Scalar | minimalBoundaryPermeability |
Minimal limit for the boundary permeability. More... | |
Scalar | ErrorTermFactor_ |
Handling of error term: relaxation factor. More... | |
Scalar | ErrorTermLowerBound_ |
Handling of error term: lower bound for error dampening. More... | |
Scalar | ErrorTermUpperBound_ |
Matrix | A_ |
Global stiffnes matrix (sparse matrix which is build by the initializeMatrix() function) More... | |
RHSVector | f_ |
Right hand side vector. More... | |
Static Protected Attributes | |
static constexpr int | pressureType = getPropValue<TypeTag, Properties::PressureFormulation>() |
gives kind of pressure used ( \( 0 = p_w \), \( 1 = p_n \), \( 2 = p_{global} \)) More... | |
general methods for output | |
TransportSolutionType | updateEstimate_ |
Update estimate for changes in volume for the pressure equation. More... | |
VtkMultiWriter< GridView > | initializationOutputWriter_ |
output for the initialization procedure More... | |
Scalar | maxError_ |
Maximum volume error of all cells. More... | |
Scalar | incp_ |
Increment for the volume derivative w.r.t pressure. More... | |
void | initializationOutput (double pseudoTS=0.) |
Write additional debug info in a special writer. More... | |
|
protected |
|
protected |
|
protected |
Indices of matrix and rhs entries.
During the assembling of the global system of equations get-functions are called (getSource(), getFlux(), etc.), which return global matrix or right hand side entries in a vector. These can be accessed using following indices:
Enumerator | |
---|---|
rhs | index for the right hand side entry |
matrix | index for the global matrix entry |
|
protectedinherited |
Indices of matrix and rhs entries.
During the assembling of the global system of equations get-functions are called (getSource(), getFlux(), etc.), which return global matrix or right hand side entries in a vector. These can be accessed using following indices:
Enumerator | |
---|---|
rhs | index for the right hand side entry |
matrix | index for the global matrix entry |
|
inline |
Constructs a FVPressure2P2CPC object.
problem | a problem class object |
|
inline |
Write data files.
writer | The writer |
void Dumux::FVPressure2P2CMultiPhysics< TypeTag >::assemble | ( | bool | first | ) |
function which assembles the system of equations to be solved
for first == true, this function assembles the matrix and right hand side for the solution of the pressure field in the same way as in the class FVPressure2P. for first == false, the approach is changed to
\[-\frac{\partial V}{\partial p} \frac{\partial p}{\partial t}+\sum_{\kappa}\frac{\partial V}{\partial m^{\kappa}}\nabla\cdot \left(\sum_{\alpha}C_{\alpha}^{\kappa}\mathbf{v}_{\alpha}\right) =\sum_{\kappa}\frac{\partial V}{\partial m^{\kappa}}q^{\kappa} \]
. See Paper SPE 99619. This is done to account for the volume effects which appear when gas and liquid are dissolved in each other.
first | Flag if pressure field is unknown |
|
inlineinherited |
|
inline |
Function for deserialization of the pressure field.
Function needed for restart option. Reads the pressure of a grid element from a restart file.
instream | Stream from the restart file. |
element | Grid element |
void Dumux::FVPressure2P2CMultiPhysics< TypeTag >::get1pFlux | ( | EntryType & | entries, |
const Intersection & | intersection, | ||
const CellData & | cellDataI | ||
) |
The compositional single-phase flux in the multiphysics framework.
If only single-phase conditions are encountered, the flux expression simplifies to (written for the case where the wetting phase is only present):
\[ A_{\gamma} \mathbf{n}_{\gamma}^T \mathbf{K} \lambda_w \mathbf{d}_{ij} \left( \frac{p_{w,j}^t - p^{t}_{w,i}}{\Delta x} + \varrho_{w} \mathbf{g}^T \mathbf{d}_{ij} \right) . \]
entries | The Matrix and RHS entries |
intersection | Intersection between cell I and J |
cellDataI | Data of cell I |
void Dumux::FVPressure2P2CMultiPhysics< TypeTag >::get1pFluxOnBoundary | ( | EntryType & | entries, |
const Intersection & | intersection, | ||
const CellData & | cellDataI | ||
) |
The compositional single-phase flux in the multiphysics framework.
If only single-phase conditions are encountered, the flux expression simplifies to (written for the case where the wetting phase is only present):
\[ A_{\gamma} \mathbf{n}_{\gamma}^T \mathbf{K} \varrho_w \lambda_w \mathbf{d}_{i-Boundary} \left( \frac{p_{w,Boundary}^t - p^{t}_{w,i}}{\Delta x} + \varrho_{w} \mathbf{g}^T \mathbf{d}_{i-Boundary} \right) . \]
If a Neumann BC is set, the given (mass-)flux is directly multiplied by the volume derivative and inserted.
entries | The Matrix and RHS entries |
intersection | Intersection between cell I and J |
cellDataI | Data of cell I |
void Dumux::FVPressure2P2CMultiPhysics< TypeTag >::get1pSource | ( | EntryType & | sourceEntry, |
const Element & | elementI, | ||
const CellData & | cellDataI | ||
) |
Assembles the source term.
The source is translated into a volumentric source term:
\[ V_i \sum_{\kappa} \frac{1}{\varrho} q^{\kappa}_i \; , \]
because under singlephase conditions
\[ \frac{\partial v_{t}}{\partial C^{\kappa}} \approx \frac{1}{\varrho} \]
.
sourceEntry | The Matrix and RHS entries |
elementI | The element I |
cellDataI | Data of cell I |
void Dumux::FVPressure2P2CMultiPhysics< TypeTag >::get1pStorage | ( | EntryType & | storageEntry, |
const Element & | elementI, | ||
CellData & | cellDataI | ||
) |
Assembles the storage term for a 1p cell in a multiphysics framework.
The storage term comprises the (single-phase) compressibility (due to a change in pressure from last timestep):
\[ V_i c_{i} \frac{p^t_i - p^{t-\Delta t}_i}{\Delta t} \]
and the damped error introduced by the incorrect transport of the last timestep:
\[ V_i \alpha_r \frac{v_{t} - \phi}{\Delta t} \]
. The latter is damped according to Fritz 2011.
storageEntry | The Matrix and RHS entries |
elementI | The element I |
cellDataI | Data of cell I |
|
inherited |
Get flux at an interface between two cells.
for first == true, the flux is calculated in traditional fractional-flow forn as in FVPressure2P. for first == false, the flux thorugh \( \gamma \) is calculated via a volume balance formulation
\[ - A_{\gamma} \mathbf{n}^T_{\gamma} \mathbf{K} \sum_{\alpha} \varrho_{\alpha} \lambda_{\alpha} \mathbf{d}_{ij} \left( \frac{p_{\alpha,j}^t - p^{t}_{\alpha,i}}{\Delta x} + \varrho_{\alpha} \mathbf{g}^T \mathbf{d}_{ij} \right) \sum_{\kappa} X^{\kappa}_{\alpha} \frac{\partial v_{t}}{\partial C^{\kappa}} + V_i \frac{A_{\gamma}}{U_i} \mathbf{d}^T \mathbf{K} \sum_{\alpha} \varrho_{\alpha} \lambda_{\alpha} \mathbf{d}_{ij} \left( \frac{p_{\alpha,j}^t - p^{t}_{\alpha,i}}{\Delta x} + \varrho_{\alpha} \mathbf{g}^T \mathbf{d}_{ij} \right) \sum_{\kappa} X^{\kappa}_{\alpha} \frac{\frac{\partial v_{t,j}}{\partial C^{\kappa}_j}-\frac{\partial v_{t,i}}{\partial C^{\kappa}_i}}{\Delta x} \]
This includes a boundary integral and a volume integral, because \( \frac{\partial v_{t,i}}{\partial C^{\kappa}_i} \) is not constant. Here, \( \mathbf{d}_{ij} \) is the normalized vector connecting the cell centers, and \( \mathbf{n}_{\gamma} \) represents the normal of the face \( \gamma \).
entries | The Matrix and RHS entries |
intersection | Intersection between cell I and J |
cellDataI | Data of cell I |
first | Flag if pressure field is unknown |
|
inherited |
Get flux on Boundary.
for first == true, the flux is calculated in traditional fractional-flow forn as in FVPressure2P. for first == false, the flux thorugh \( \gamma \) is calculated via a volume balance formulation
\[ - A_{\gamma} \mathbf{n}^T_{\gamma} \mathbf{K} \sum_{\alpha} \varrho_{\alpha} \lambda_{\alpha} \mathbf{d}_{ij} \left( \frac{p_{\alpha,j}^t - p^{t}_{\alpha,i}}{\Delta x} + \varrho_{\alpha} \mathbf{g}^T \mathbf{d}_{ij} \right) \sum_{\kappa} \frac{\partial v_{t}}{\partial C^{\kappa}} X^{\kappa}_{\alpha} \;, \]
where we skip the volume integral assuming \( \frac{\partial v_{t,i}}{\partial C^{\kappa}_i} \) to be constant at the boundary. Here, \( \mathbf{d}_{ij} \) is the normalized vector connecting the cell centers, and \( \mathbf{n}_{\gamma} \) represents the normal of the face \( \gamma \).
If a Neumann BC is set, the given (mass-)flux is directly multiplied by the volume derivative and inserted.
entries | The Matrix and RHS entries |
intersection | Intersection between cell I and J |
cellDataI | Data of cell I |
first | Flag if pressure field is unknown |
|
inherited |
Assembles the source term.
for first == true, a source is implemented as in FVPressure2P. for first == false, the source is translated into a volumentric source term:
\[ V_i \sum_{\kappa} \frac{\partial v_{t}}{\partial C^{\kappa}} q^{\kappa}_i \]
.
sourceEntry | The Matrix and RHS entries |
elementI | The element I |
cellDataI | Data of cell I |
first | Flag if pressure field is unknown |
|
inherited |
Assembles the storage term.
for first == true, there is no storage contribution. for first == false, the storage term comprises the compressibility (due to a change in pressure from last timestep):
\[ V_i c_{t,i} \frac{p^t_i - p^{t-\Delta t}_i}{\Delta t} \]
and the damped error introduced by the incorrect transport of the last timestep:
\[ V_i \alpha_r \frac{v_{t} - \phi}{\Delta t} \]
. The latter is damped according to Fritz 2011.
storageEntry | The Matrix and RHS entries |
elementI | The element I |
cellDataI | Data of cell I |
first | Flag if pressure field is unknown |
|
inlineinherited |
Returns the global matrix of the last pressure solution step.
|
inlineinherited |
Write additional debug info in a special writer.
To visualize the different steps through the initialization procedure, we use very small pseudo time steps only for the writer! This is only for debugging of the initialization procedure.
pseudoTS | Time steps that only appear in the writer, not real. |
|
inlineinherited |
Initialize pressure model.
Function initializes the sparse matrix to solve the global system of equations and sets/calculates the initial pressure
|
inline |
|
protectedinherited |
Initialize the global matrix of the system of equations to solve.
|
protectedinherited |
Initialize the global matrix of the system of equations to solve.
|
protectedinherited |
Initialize the global matrix of the system of equations to solve.
|
inlineprotectedinherited |
Initialization of the pressure solution vector: Initialization with meaningful values may.
|
inherited |
initializes the fluid distribution and hereby the variables container
It differs from updateMaterialLaws() because there are two possible initial conditions: saturations and concentration.
compositional | flag that determines if compositional effects are regarded, i.e. a reasonable pressure field is known with which compositions can be calculated. |
|
inlineprotectedinherited |
Returns the vector containing the pressure solution.
|
inlineprotectedinherited |
Returns the vector containing the pressure solution.
|
inlineinherited |
Public access function for the primary pressure variable.
Function returns the cell pressure value at index eIdxGlobal
eIdxGlobal | Global index of a grid cell |
|
inlineinherited |
|
inlineinherited |
Returns the right hand side vector of the last pressure solution step.
|
inline |
Function for serialization of the pressure field.
Function needed for restart option. Writes the pressure of a grid element to a restart file.
outstream | Stream into the restart file. |
element | Grid element |
|
inlineinherited |
Set a pressure to be fixed at a certain cell.
Allows to fix a pressure somewhere (at one certain cell) in the domain. This can be necessary e.g. if only Neumann boundary conditions are defined. The pressure is fixed until the unsetFixPressureAtIndex()
function is called
pressure | Pressure value at eIdxGlobal |
eIdxGlobal | Global index of a grid cell |
|
protectedinherited |
Solves the global system of equations to get the spatial distribution of the pressure.
|
inlineinherited |
Reset the fixed pressure state.
No pressure is fixed inside the domain until setFixPressureAtIndex()
function is called again.
eIdxGlobal | Global index of a grid cell |
|
inlineinherited |
Compositional pressure solution routine: update estimate for secants, assemble, solve.
An update estime (transport step acoording to old pressure field) determines changes in mass, composition, which is used to calculate volume derivatives entering the pressure equation, as well as an approximate guess for time step size for the storage terms in the p.e. Afterwards, the system is assembled and solved for pressure.
void Dumux::FVPressure2P2CMultiPhysics< TypeTag >::update1pMaterialLawsInElement | ( | const Element & | elementI, |
CellData & | cellData, | ||
bool | postTimeStep | ||
) |
updates secondary variables of one single phase cell
For each element, the secondary variables are updated according to the primary variables. Only a simple flash calulation has to be carried out, as phase distribution is already known: single-phase.
elementI | The element |
cellData | The cell data of the current element |
postTimeStep | Flag indicating if we have just completed a time step |
void Dumux::FVPressure2P2CMultiPhysics< TypeTag >::updateMaterialLaws | ( | bool | postTimeStep = false | ) |
constitutive functions are updated once if new concentrations are calculated and stored in the variables container
In contrast to the standard sequential 2p2c model ( FVPressure2P2C<TypeTag>::updateMaterialLaws() ), this method also holds routines to adapt the subdomain. The subdomain indicates weather we are in 1p domain (value = 1) or in the two phase subdomain (value = 2). Note that the type of flash, i.e. the type of FluidState (FS), present in each cell does not have to coincide with the subdomain. If a cell will be simple and was complex, a complex FS is available, so next time step will use this complex FS, but updateMaterialLaw afterwards will finally transform that to simple FS.
postTimeStep | Flag indicating method is called from Problem::postTimeStep() |
|
inherited |
Updates secondary variables of one cell.
For each element, the secondary variables are updated according to the primary variables. In case the method is called after the Transport, i.e. at the end / post time step, CellData2p2c.reset() resets the volume derivatives for the next time step.
element | The element |
postTimeStep | Flag indicating if we have just completed a time step |
|
inlineinherited |
|
inherited |
Partial derivatives of the volumes w.r.t. changes in total concentration and pressure.
This method calculates the volume derivatives via a secant method, where the secants are gained in a pre-computational step via the transport equation and the last TS size. The partial derivatives w.r.t. mass are defined as \( \frac{\partial v}{\partial C^{\kappa}} = \frac{\partial V}{\partial m^{\kappa}}\)
globalPos | The global position of the current element |
element | The current element |
|
protectedinherited |
Global stiffnes matrix (sparse matrix which is build by the initializeMatrix()
function)
|
protectedinherited |
Enables the volume integral of the pressure equation.
|
protectedinherited |
Handling of error term: relaxation factor.
|
protectedinherited |
Handling of error term: lower bound for error dampening.
|
protectedinherited |
Handling of error term: upper bound for error dampening
|
protectedinherited |
Right hand side vector.
|
protected |
vector including the gravity constant
|
protectedinherited |
Increment for the volume derivative w.r.t pressure.
|
protectedinherited |
output for the initialization procedure
|
protectedinherited |
Maximum volume error of all cells.
|
protectedinherited |
Minimal limit for the boundary permeability.
|
protected |
vector holding next subdomain
|
staticconstexprprotected |
gives kind of pressure used ( \( 0 = p_w \), \( 1 = p_n \), \( 2 = p_{global} \))
|
protectedinherited |
|
protectedinherited |
Enables regulation of permeability in the direction of a Dirichlet Boundary Condition.
|
protected |
A timer for the time spent on the multiphysics framework.
|
protectedinherited |
Update estimate for changes in volume for the pressure equation.