Loading...
Searching...
No Matches
Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct > Class Template Reference

Solve mixed formulation of incompressible flow modelled by Darcy's law. More...

`#include <IncompFlowSolverHybrid.hpp>`

## Public Types

typedef const FlowSolution & | SolutionType
| Type representing the solution to the problem defined by the parameters to.

## Public Member Functions

template<class Point>
void | init (const GridInterface &g, const RockInterface &r, const Point &grav, const BCInterface &bc)
| All-in-one initialization routine.
void | clear ()
| Clear all topologic, geometric and rock-dependent information currently held in internal data structures.
void | initSystemStructure (const GridInterface &g, const BCInterface &bc)
| Compute structure of coefficient matrix in final system of linear equations for this flow problem.
template<class Point>
void | computeInnerProducts (const RockInterface &r, const Point &grav)
| Compute static (i.e., independent of saturation) parts of the spatially varying inner products for each grid cell.
template<class FluidInterface>
void | solve (const FluidInterface &r, const std::vector< double > &sat, const BCInterface &bc, const std::vector< double > &src, double residual_tolerance=1e-8, int linsolver_verbosity=1, int linsolver_type=1, bool same_matrix=false, int linsolver_maxit=0, double prolongate_factor=1.6, int smooth_steps=1)
| Construct and solve system of linear equations for the pressure values on each interface/contact between neighbouring grid cells.
double | postProcessFluxes ()
| Postprocess the solution fluxes.
SolutionType | getSolution ()
| Recover the solution to the problem defined by the parameters to method.
template<typename charT, class traits>
void | printStats (std::basic_ostream< charT, traits > &os)
| Print statistics about the connections in the current model.
void | printSystem (const std::string &prefix)
| Output current system of linear equations to permanent storage in files.

## Detailed Description

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>
class Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >

Solve mixed formulation of incompressible flow modelled by Darcy's law.

@f] The solver is based on a hybrid discretization scheme which yields a system of linear equations of the form

=

@f] where represents the interface fluxes for each interface for each cell, are the cell pressures and are the interface pressures.

Through a Schur complement analysis, the above system of linear equations is reduced to an equivalent system of the form

=

@f] where , , and . Similarly, and .

The system is the system of linear equations which is actually solved using separte linear system solver software. Then, a back substitution process yields the cell pressures and interface fluxes by solving the simpler systems

@f] Specifically, the matrix is diagonal (a single scalar per grid cell) and the matrix is available as a bi-product of the Schur reduction process. Finally, in the case of mimetic discretizations using the

MimeticIPEvaluator
Opm::MimeticIPEvaluator
Definition MimeticIPEvaluator.hpp:85

class, the matrix is directly available through an explicit formula.

The matrix is a product of two parts–one which depends only on the grid (i.e., the topology and geomtry) and the geological properties (i.e., the permeability field) and one part which depends only on the temporally varying fluid data (i.e., viscosities and saturation dependent relative permeabilities). Furthermore, the matrices and depend only on the grid topolgy. Consequently, the cost of assembling the system may be reduced, at least for simulations involving multiple pressure solves, by performing all of the static (i.e., the parts dependent only on the grid and geological properties) assembly once, at system initalization. Then the system assembly process consists of calculating the saturation dependent contributions and adding these into the system of linear equations.

The code is structured in a manner similar to traditional FEM software implementations. In particular, we assemble the coefficient matrix and system right hand side on a cell-by-cell basis. This yields a number of important simplifications. In particular, the matrix reduces to the identity on a single cell. Its effects are produced only through a local-to-global map of degrees of freedom. Similarly, on a single grid cell, the matrix with the number of elements equal to the number of interfaces (i.e., neighbours) of the cell. Finally, the matrix is computed as

where denotes the static part of and denotes the total fluid mobility in the grid cell. For ease of implementation of the back substitution process, we compute and store the values of , and during each system assembly process.

A final quirk of the implementation is that we always solve for every interface pressure, even if a given pressure value is known through, e.g., a prescribed pressure boundary condition. This feature enables changing the type and value of a set of boundary conditions between pressure solves without having to reconstruct the sparsity structure of the coefficient matrix –at least while no connections are introduced or removed. Periodic boundary condtions introduce new connections (i.e., new non-zero coefficient matrix entries) so in order to gain maximum efficiency, you should initialize the flow solver system only after you know the existence of all periodic boundary conditions which will be employed in a given simulation.

The use of non-periodic boundary conditions after periodic connections have been introduced is fully supported. This mode introduces a performance caveat. Specifically, a few entries which would otherwise not have been entered into the matrix due to being automatically zero will become explicit zeros instead. Consequently, each iteration of an iterative solver for the system becomes slightly more expensive as explicit zero-multiplications occur. Moreover, introducing explicit non-zero representation of matrix entries which are always zero increases the demand for memory resources.

The flow solver is initialized in a three-step process. The first step, represented by method

clear()
Opm::IncompFlowSolverHybrid::clear
void clear()
Clear all topologic, geometric and rock-dependent information currently held in internal data structu...
Definition IncompFlowSolverHybrid.hpp:497

, releases any previously held data in the internal data structures and prepares the solver system for defining a new problem. The second step, represented by the

initSystemStructure()
Opm::IncompFlowSolverHybrid::initSystemStructure
void initSystemStructure(const GridInterface &g, const BCInterface &bc)
Compute structure of coefficient matrix in final system of linear equations for this flow problem.
Definition IncompFlowSolverHybrid.hpp:533

enumerates the primary degrees of freedom (i.e., the interface pressure values) and determines, and allocates, the coefficient matrix non-zero structure. The third and final step, represented by method

computeInnerProducts()
Opm::IncompFlowSolverHybrid::computeInnerProducts
void computeInnerProducts(const RockInterface &r, const Point &grav)
Compute static (i.e., independent of saturation) parts of the spatially varying inner products for ...
Definition IncompFlowSolverHybrid.hpp:565

, computes the static (i.e., geology and geomtry-dependent) inner product matrices . Method

computeInnerProducts()

is offered separately in order to support solve several different property models on the same underlying geometry (grid). Finally, method

init()
Opm::IncompFlowSolverHybrid::init
void init(const GridInterface &g, const RockInterface &r, const Point &grav, const BCInterface &bc)
All-in-one initialization routine.
Definition IncompFlowSolverHybrid.hpp:476

is a convenience method which calls the other initialization methods in sequence.

Following solver intialization, a sequence of flow problems differing only by boundary condition type/value and/or differing fluid saturation values may be resolved by separate calls to the

solve()
Opm::IncompFlowSolverHybrid::solve
void solve(const FluidInterface &r, const std::vector< double > &sat, const BCInterface &bc, const std::vector< double > &src, double residual_tolerance=1e-8, int linsolver_verbosity=1, int linsolver_type=1, bool same_matrix=false, int linsolver_maxit=0, double prolongate_factor=1.6, int smooth_steps=1)
Construct and solve system of linear equations for the pressure values on each interface/contact betw...
Definition IncompFlowSolverHybrid.hpp:658

method. At any time following a call to

solve()

may the solution, represented by the type

SolutionType
Opm::IncompFlowSolverHybrid::SolutionType
const FlowSolution & SolutionType
Type representing the solution to the problem defined by the parameters to.
Definition IncompFlowSolverHybrid.hpp:811

, to the most recently resolved flow problem be retrieved through the

getSolution()
Opm::IncompFlowSolverHybrid::getSolution
SolutionType getSolution()
Recover the solution to the problem defined by the parameters to method.
Definition IncompFlowSolverHybrid.hpp:821

method. We note that

SolutionType

is always a reference-to-const.

Template Parameters

  GridInterface

  Type presenting an interface to a grid (typically a discretized geological model). The type is assumed to expose a forward iterator type,
  CellIter
  , and a pair of delimiters
  cellbegin()
  and
  cellend()
  to traverse the cells of a grid. The cell iterator, in turn, is expected to expose a forward iterator type
  FaceIter
  , and a pair of delimiters
  facebegin()
  and
  faceend()
  to traverse the faces of a cell.

  RockInterface | Type presenting an interface to reservoir properties such as permeability and porosity. The type is assumed to expose a method, permeability() through which the assigned permeability of a single grid cell, represented by a GridInterface::CellIter , may be recovered.

  BCInterface

  Type presenting an interface to boundary conditions. The type is expected to provide a method,
  flowCond()
  , from which a
  BCInterface::FlowBC
  boundary condtion type of any given boundary face may be recovered. The flow boundary condition type is expected to provide methods for quering the type of boundary condition (i.e., prescribed pressure values, prescribed flux values, or prescribed pressure drops in the case of periodic boundary conditions) as well as the numerical values of these conditions.

  InnerProduct | Type presenting a specific inner product defining a discretization of the Darcy equation.

## Member Typedef Documentation

## SolutionType

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>

typedef const FlowSolution& Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::SolutionType

Type representing the solution to the problem defined by the parameters to.

solve()

. Always a reference-to-const. The

SolutionType

exposes methods

pressure(c)

and

outflux(f)

from which the cell pressure in cell

*c

and outward-pointing flux across interface

*f

may be recovered.

## Member Function Documentation

## clear()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>

void Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::clear | ( |  | ) |

inline

Clear all topologic, geometric and rock-dependent information currently held in internal data structures.

Release all memory resources in preparation of constructing a solver for a new problem. Method

clear()

must be called prior to any other method of the class.

## computeInnerProducts()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>
template<class Point>

void Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::computeInnerProducts | ( | const RockInterface & | r,
|  | const Point & | grav )

inline

Compute static (i.e., independent of saturation) parts of the spatially varying inner products for each grid cell.

This is often a fairly expensive process, so in order to increase the speed of the system assembly we do not wish to perform this computation for each time step unless we have to. Moreover, if the permeability field changes but the grid connections remain static, we do not have to build a new matrix structure but need only compute new values for these static inner products.

Parameters

  [in] | r | The reservoir properties of each grid cell. In method computeInnerProducts() , we only inspect the permeability field of r .
  [in] | grav | Gravity

## getSolution()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>

SolutionType Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::getSolution | ( |  | ) |

inline

Recover the solution to the problem defined by the parameters to method.

solve()

. This solution is meaningless without a previous call to method

solve()

.

Returns
  The current solution.

## init()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>
template<class Point>

void Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::init | ( | const GridInterface & | g,
|  | const RockInterface & | r,
|  | const Point & | grav,
|  | const BCInterface & | bc )

inline

All-in-one initialization routine.

Enumerates all grid connections, allocates sufficient space, defines the structure of the global system of linear equations for the contact pressures, and computes the permeability dependent inner products for all of the grid's cells.

Parameters

  [in] | g | The grid.
  [in] | r | The reservoir properties of each grid cell.
  [in] | grav | Gravity vector. Its Euclidian two-norm value represents the strength of the gravity field (in units of m/s^2) while its direction is the direction of gravity in the current model.
  [in] | bc | The boundary conditions describing how the current flow problem interacts with the outside world. This is used only for the purpose of introducing additional couplings in the case of periodic boundary conditions. The specific values of the boundary conditions are not inspected in init() .

## initSystemStructure()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>

void Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::initSystemStructure | ( | const GridInterface & | g,
|  | const BCInterface & | bc )

inline

Compute structure of coefficient matrix in final system of linear equations for this flow problem.

Allocates all storage needed by the flow solver itself.

Parameters

  [in] | g | The grid.
  [in] | bc | The boundary conditions describing how the current flow problem interacts with the outside world. This is used only for the purpose of introducing additional couplings in the case of periodic boundary conditions. The specific values of the boundary conditions are not inspected in init() .

## postProcessFluxes()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>

double Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::postProcessFluxes | ( |  | ) |

inline

Postprocess the solution fluxes.

This method modifies the solution object so that out-fluxes of twin faces (that is, the two faces on a cell-cell intersection) will be made antisymmetric.

Returns
  The maximum modification made to the fluxes.

## printStats()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>
template<typename charT, class traits>

void Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::printStats | ( | std::basic_ostream< charT, traits > & | os | ) |

inline

Print statistics about the connections in the current model.

This is mostly for debugging purposes and should only rarely be used in client code.

Template Parameters

  charT | Character type of output stream.
  traits | Character traits of charT .

Parameters

  os | Output stream into which the statistics will be printed.

## printSystem()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>

void Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::printSystem | ( | const std::string & | prefix | ) |

inline

Output current system of linear equations to permanent storage in files.

One file for the coefficient matrix and one file for the right hand side. This is mostly useful whilst debugging the solver and should rarely be used from client code.

The system is stored in a format which is suitable for importing into MATLAB or MATLAB-compatible software. In particular, using the MATLAB

LOAD

and

SPCONVERT

functions makes it easy to reconstruct the system of linear equations from within MATLAB.

Parameters

  [in] | prefix | Prefix from which file names for the coefficient matrix and right hand side data are derived. Specifically, the matrix data is output to the file prefix + "-mat.dat" while the right hand side data is output to the file prefix + "-rhs.dat" .

## solve()

template<class GridInterface, class RockInterface, class BCInterface, template< class GridIF, class RockIF > class InnerProduct>
template<class FluidInterface>

void Opm::IncompFlowSolverHybrid< GridInterface, RockInterface, BCInterface, InnerProduct >::solve | ( | const FluidInterface & | r,
|  | const std::vector< double > & | sat,
|  | const BCInterface & | bc,
|  | const std::vector< double > & | src,
|  | double | residual_tolerance = 1e-8,
|  | int | linsolver_verbosity = 1,
|  | int | linsolver_type = 1,
|  | bool | same_matrix = false,
|  | int | linsolver_maxit = 0,
|  | double | prolongate_factor = 1.6,
|  | int | smooth_steps = 1 )

inline

Construct and solve system of linear equations for the pressure values on each interface/contact between neighbouring grid cells.

Recover cell pressure and interface fluxes. Following a call to

solve()
@encode, you may recover the flow solution from the
@code getSolution()

method.

Template Parameters

  FluidInterface | Type presenting an interface to fluid properties such as density, mobility &c. The type is expected to provide methods phaseMobilities() and phaseDensities() for phase mobility and density in a single cell, respectively.

Parameters

  [in] | r | The fluid properties of each grid cell. In method solve() we query this object for the phase mobilities (i.e., r.phaseMobilities() ) and the phase densities (i.e., phaseDensities() ) of each phase.
  [in] | sat | Saturation of primary phase. One scalar value for each grid cell. This parameter currently limits IncompFlowSolverHybrid Opm::IncompFlowSolverHybrid Solve mixed formulation of incompressible flow modelled by Darcy's law. Definition IncompFlowSolverHybrid.hpp:365 to two-phase flow problems.
  [in] | bc | The boundary conditions describing how the current flow problem interacts with the outside world. Method solve() inspects the actual values of the boundary conditions whilst forming the system of linear equations.

Specifically, the

bc.flowCond(bid)

method is expected to yield a valid

FlowBC
Opm::FlowBC
A class for representing a flow boundary condition.
Definition BoundaryConditions.hpp:122

object for which the methods

pressure()

,

pressureDifference()

, and

outflux()

yield valid responses depending on the type of the object.

Parameters

  [in] | src | Explicit source terms. One scalar value for each grid cell representing the rate (in units of m^3/s) of fluid being injected into (>0) or extracted from (<0) a given grid cell.
  [in] | residual_tolerance | Control parameter for iterative linear solver software. The iteration process is terminated when the norm of the linear system residual is less than residual_tolerance times the initial residual.
  [in] | linsolver_verbosity | Control parameter for iterative linear solver software. Verbosity level 0 prints nothing, level 1 prints summary information, level 2 prints data for each iteration.
  [in] | linsolver_type | Control parameter for iterative linear solver software. Type 0 selects a ILU0/CG solver, type 1 selects AMG/CG, type 2 selects KAMG/CG, type 3 selects AMG/CG with fast Gauss-Seidel smoothing
  [in] | linsolver_maxit | maximum iterations allowed
  [in] | prolongate_factor | Factor to scale the prolongated coarse coarse grid correction
  [in] | smooth_steps | Number of smoothing steps to be used for pre and post smoothing in AMG
  [in] | same_matrix | Whether the matrix is the same as in the previous solve.

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

- opm/porsol/mimetic/IncompFlowSolverHybrid.hpp