My Project
Loading...
Searching...
No Matches
Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar > Class Template Reference

Public Types

enum class  TransUpdateQuantities { Trans , All }
 Compute all transmissibilities. More...
 
using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>
 
using DimVector = Dune::FieldVector<Scalar, dimWorld>
 

Public Member Functions

 Transmissibility (const EclipseState &eclState, const GridView &gridView, const CartesianIndexMapper &cartMapper, const Grid &grid, std::function< std::array< double, dimWorld >(int)> centroids, bool enableEnergy, bool enableDiffusivity, bool enableDispersivity)
 
const DimMatrix & permeability (unsigned elemIdx) const
 Return the permeability for an element.
 
Scalar transmissibility (unsigned elemIdx1, unsigned elemIdx2) const
 Return the transmissibility for the intersection between two elements.
 
Scalar transmissibilityBoundary (unsigned elemIdx, unsigned boundaryFaceIdx) const
 Return the transmissibility for a given boundary segment.
 
Scalar thermalHalfTrans (unsigned insideElemIdx, unsigned outsideElemIdx) const
 Return the thermal "half transmissibility" for the intersection between two elements.
 
Scalar thermalHalfTransBoundary (unsigned insideElemIdx, unsigned boundaryFaceIdx) const
 
Scalar diffusivity (unsigned elemIdx1, unsigned elemIdx2) const
 Return the diffusivity for the intersection between two elements.
 
Scalar dispersivity (unsigned elemIdx1, unsigned elemIdx2) const
 Return the dispersivity for the intersection between two elements.
 
void finishInit (const std::function< unsigned int(unsigned int)> &map={})
 Actually compute the transmissibility over a face as a pre-compute step.
 
void update (bool global, TransUpdateQuantities update_quantities=TransUpdateQuantities::All, const std::function< unsigned int(unsigned int)> &map={}, bool applyNncMultRegT=false)
 

Protected Member Functions

void updateFromEclState_ (bool global)
 
void removeNonCartesianTransmissibilities_ (bool removeAll)
 
void applyAllZMultipliers_ (Scalar &trans, unsigned insideFaceIdx, unsigned outsideFaceIdx, unsigned insideCartElemIdx, unsigned outsideCartElemIdx, const TransMult &transMult, const std::array< int, dimWorld > &cartDims)
 Apply the Multipliers for the case PINCH(4)==TOPBOT.
 
std::array< std::vector< double >, 3 > createTransmissibilityArrays_ (const std::array< bool, 3 > &is_tran)
 Creates TRANS{XYZ} arrays for modification by FieldProps data.
 
void resetTransmissibilityFromArrays_ (const std::array< bool, 3 > &is_tran, const std::array< std::vector< double >, 3 > &trans)
 overwrites calculated transmissibilities
 
template<class Intersection >
void computeFaceProperties (const Intersection &intersection, const int, const int, const int, const int, DimVector &faceCenterInside, DimVector &faceCenterOutside, DimVector &faceAreaNormal, std::false_type) const
 
template<class Intersection >
void computeFaceProperties (const Intersection &intersection, const int insideElemIdx, const int insideFaceIdx, const int outsideElemIdx, const int outsideFaceIdx, DimVector &faceCenterInside, DimVector &faceCenterOutside, DimVector &faceAreaNormal, std::true_type) const
 
void applyNncToGridTrans_ (const std::unordered_map< std::size_t, int > &cartesianToCompressed)
 
void applyPinchNncToGridTrans_ (const std::unordered_map< std::size_t, int > &cartesianToCompressed)
 Applies the previous calculate transmissibilities to the NNCs created via PINCH.
 
void applyEditNncToGridTrans_ (const std::unordered_map< std::size_t, int > &globalToLocal)
 Multiplies the grid transmissibilities according to EDITNNC.
 
void applyEditNncrToGridTrans_ (const std::unordered_map< std::size_t, int > &globalToLocal)
 Resets the grid transmissibilities according to EDITNNCR.
 
void applyNncMultreg_ (const std::unordered_map< std::size_t, int > &globalToLocal)
 
void applyEditNncToGridTransHelper_ (const std::unordered_map< std::size_t, int > &globalToLocal, const std::string &keyword, const std::vector< NNCdata > &nncs, const std::function< KeywordLocation(const NNCdata &)> &getLocation, const std::function< void(Scalar &, const Scalar &)> &apply)
 
void extractPermeability_ ()
 
void extractPermeability_ (const std::function< unsigned int(unsigned int)> &map)
 
void extractPorosity_ ()
 
void extractDispersion_ ()
 
void computeHalfTrans_ (Scalar &halfTrans, const DimVector &areaNormal, int faceIdx, const DimVector &distance, const DimMatrix &perm) const
 
void computeHalfDiffusivity_ (Scalar &halfDiff, const DimVector &areaNormal, const DimVector &distance, const Scalar &poro) const
 
DimVector distanceVector_ (const DimVector &faceCenter, const unsigned &cellIdx) const
 
void applyMultipliers_ (Scalar &trans, unsigned faceIdx, unsigned cartElemIdx, const TransMult &transMult) const
 
void applyNtg_ (Scalar &trans, unsigned faceIdx, unsigned elemIdx, const std::vector< double > &ntg) const
 

Protected Attributes

std::vector< DimMatrix > permeability_
 
std::vector< Scalar > porosity_
 
std::vector< Scalar > dispersion_
 
std::unordered_map< std::uint64_t, Scalar > trans_
 
const EclipseState & eclState_
 
const GridView & gridView_
 
const CartesianIndexMapper & cartMapper_
 
const Grid & grid_
 
std::function< std::array< double, dimWorld >(int)> centroids_
 
Scalar transmissibilityThreshold_
 
std::map< std::pair< unsigned, unsigned >, Scalar > transBoundary_
 
std::map< std::pair< unsigned, unsigned >, Scalar > thermalHalfTransBoundary_
 
bool enableEnergy_
 
bool enableDiffusivity_
 
bool enableDispersivity_
 
bool warnEditNNC_ = true
 
std::unordered_map< std::uint64_t, Scalar > thermalHalfTrans_
 
std::unordered_map< std::uint64_t, Scalar > diffusivity_
 
std::unordered_map< std::uint64_t, Scalar > dispersivity_
 
const LookUpData< Grid, GridView > lookUpData_
 
const LookUpCartesianData< Grid, GridView > lookUpCartesianData_
 

Member Enumeration Documentation

◆ TransUpdateQuantities

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
enum class Opm::Transmissibility::TransUpdateQuantities
strong

Compute all transmissibilities.

Parameters
[in]globalWhether or not to call update() on all processes. Also, this updates the "thermal half transmissibilities" if energy is enabled.
[in]transIndicating whether we only allocate and upate trans_ without considering thermalHalfTrans_, diffusivity_, dispersivity_. For many usage, we only need trans_, e.g. weights for domain decomposition, INIT file output. It might change following further development. Trans only update the trans_, which is related to permeability All upate rans_, thermalHalfTrans_, diffusivity_ and dispersivity_.
[in]mapUndocumented.
[in]applyNncMultRegTWhether or not to apply regional multipliers to explicit NNCs. Explicit NNCs are those entered directly in the input data, e.g., through the NNC/EDITNNC/EDITNNCR keywords, or the result of generating connections to or within numerical aquifers. Default value: false, meaning do not apply regional multipliers to explicit NNCs.

Member Function Documentation

◆ applyAllZMultipliers_()

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
void Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar >::applyAllZMultipliers_ ( Scalar & trans,
unsigned insideFaceIdx,
unsigned outsideFaceIdx,
unsigned insideCartElemIdx,
unsigned outsideCartElemIdx,
const TransMult & transMult,
const std::array< int, dimWorld > & cartDims )
protected

Apply the Multipliers for the case PINCH(4)==TOPBOT.

Parameters
pinchTopWhether PINCH(5) is TOP, otherwise ALL is assumed.

◆ applyPinchNncToGridTrans_()

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
void Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar >::applyPinchNncToGridTrans_ ( const std::unordered_map< std::size_t, int > & cartesianToCompressed)
protected

Applies the previous calculate transmissibilities to the NNCs created via PINCH.

Parameters
cartesianToCompressedVector containing the compressed index (or -1 for inactive cells) as the element at the cartesian index.

◆ computeFaceProperties()

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
template<class Intersection >
void Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar >::computeFaceProperties ( const Intersection & intersection,
const int insideElemIdx,
const int insideFaceIdx,
const int outsideElemIdx,
const int outsideFaceIdx,
DimVector & faceCenterInside,
DimVector & faceCenterOutside,
DimVector & faceAreaNormal,
std::true_type  ) const
protected

Alternatively, the actual area of the refined face can be computed as follows:

◆ createTransmissibilityArrays_()

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
std::array< std::vector< double >, 3 > Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar >::createTransmissibilityArrays_ ( const std::array< bool, 3 > & is_tran)
protected

Creates TRANS{XYZ} arrays for modification by FieldProps data.

Parameters
is_tranWhether TRAN{XYZ} will be modified. If entry is false the array will be empty
Returns
an array of vector (TRANX, TRANY, TRANZ}

◆ finishInit()

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
void Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar >::finishInit ( const std::function< unsigned int(unsigned int)> & map = {})
inline

Actually compute the transmissibility over a face as a pre-compute step.

This code actually uses the direction specific "centroids" of each element. These "centroids" are not the identical barycenter of the element, but the middle of the centers of the faces of the logical Cartesian cells, i.e., the centers of the faces of the reference elements. We do things this way because the barycenter of the element can be located outside of the element for sufficiently "ugly" (i.e., thin and "non-flat") elements which in turn leads to quite wrong permeabilities. This approach is probably not always correct either but at least it seems to be much better.

◆ resetTransmissibilityFromArrays_()

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
void Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar >::resetTransmissibilityFromArrays_ ( const std::array< bool, 3 > & is_tran,
const std::array< std::vector< double >, 3 > & trans )
protected

overwrites calculated transmissibilities

Parameters
is_tranWhether TRAN{XYZ} have been modified.
transArrays with modified transmissibilities TRAN{XYZ}

◆ thermalHalfTrans()

template<class Grid , class GridView , class ElementMapper , class CartesianIndexMapper , class Scalar >
Scalar Opm::Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar >::thermalHalfTrans ( unsigned insideElemIdx,
unsigned outsideElemIdx ) const

Return the thermal "half transmissibility" for the intersection between two elements.

The "half transmissibility" features all sub-expressions of the "thermal transmissibility" which can be precomputed, i.e. they are not dependent on the current solution:

H_t = A * (n*d)/(d*d);

where A is the area of the intersection between the inside and outside elements, n is the outer unit normal, and d is the distance between the center of the inside cell and the center of the intersection.


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