OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #3565 from blattms/centroids-vanguard

Browse Source 
[refactor] Move cell centroid lookup magic to vanguard for reuse.
This commit is contained in:
Tor Harald Sandve
2021-10-05 11:38:19 +02:00
committed by GitHub
parent 61ef539bf5 457d270bdf
commit 033995ab8f
7 changed files with 83 additions and 39 deletions
Download Patch File Download Diff File Expand all files Collapse all files

18
ebos/eclalugridvanguard.hh
View File

@@ -102,7 +102,6 @@ public:
~EclAluGridVanguard()
{
delete cartesianIndexMapper_;
delete equilCartesianIndexMapper_;
delete grid_;
delete equilGrid_;
@@ -178,6 +177,18 @@ public:
const EquilCartesianIndexMapper& equilCartesianIndexMapper() const
{ return *equilCartesianIndexMapper_; }
/*!
* \brief Get function to query cell centroids for a distributed grid.
*
* Currently this only non-empty for a loadbalanced CpGrid.
* It is a function return the centroid for the given element
* index.
*/
std::function<std::array<double,dimensionworld>(int)>
cellCentroids() const
{
return this->cellCentroids_(cartesianIndexMapper_.get());
}
protected:
void createGrids_()
{
@@ -210,7 +221,8 @@ protected:
grid_ = factory.convert(*equilGrid_, cartesianCellId_);
cartesianIndexMapper_ =
new CartesianIndexMapper(*grid_, cartesianDimension_, cartesianCellId_);
std::make_unique<CartesianIndexMapper>(*grid_, cartesianDimension_,
cartesianCellId_);
}
void filterConnections_()
@@ -222,7 +234,7 @@ protected:
EquilGrid* equilGrid_;
std::vector<int> cartesianCellId_;
std::array<int,dimension> cartesianDimension_;
CartesianIndexMapper* cartesianIndexMapper_;
std::unique_ptr<CartesianIndexMapper> cartesianIndexMapper_;
EquilCartesianIndexMapper* equilCartesianIndexMapper_;
};

42
ebos/eclbasevanguard.hh
View File

@@ -36,6 +36,7 @@
#include <opm/grid/common/GridEnums.hpp>
#include <opm/grid/common/CartesianIndexMapper.hpp>
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/NumericalAquiferCell.hpp>
#include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
#include <opm/simulators/flow/BlackoilModelParametersEbos.hpp>
#include <array>
@@ -185,6 +186,7 @@ public:
protected:
static const int dimension = Grid::dimension;
static const int dimensionworld = Grid::dimensionworld;
using Element = typename GridView::template Codim<0>::Entity;
using CartesianIndexMapper = Dune::CartesianIndexMapper<Grid>;
@@ -329,16 +331,6 @@ public:
{ return asImp_().equilCartesianIndexMapper().cartesianCoordinate(cellIdx, ijk); }
/*!
* \brief Get the cell centroids for a distributed grid.
*
* Currently this only non-empty for a loadbalanced CpGrid.
*/
const std::vector<double>& cellCentroids() const
{
return centroids_;
}
/*!
* \brief Returns the depth of a degree of freedom [m]
*
@@ -389,6 +381,36 @@ public:
}
protected:
/*!
* \brief Get function to query cell centroids for a distributed grid.
*
* Currently this only non-empty for a loadbalanced CpGrid.
* It is a function return the centroid for the given element
* index.
* \param cartMapper The cartesian index mapper for lookup of
* cartesian indices
*/
template<class CartMapper>
std::function<std::array<double,dimensionworld>(int)>
cellCentroids_(const CartMapper& cartMapper) const
{
return [this, cartMapper](int elemIdx) {
const auto& centroids = this->centroids_;
auto rank = this->gridView().comm().rank();
std::array<double,dimensionworld> centroid;
if (rank == 0) {
unsigned cartesianCellIdx = cartMapper.cartesianIndex(elemIdx);
centroid = this->eclState().getInputGrid().getCellCenter(cartesianCellIdx);
} else
{
std::copy(centroids.begin() + elemIdx * dimensionworld,
centroids.begin() + (elemIdx + 1) * dimensionworld,
centroid.begin());
}
return centroid;
};
}
void callImplementationInit()
{
asImp_().createGrids_();

14
ebos/eclcpgridvanguard.hh
View File

@@ -90,6 +90,7 @@ public:
using EquilGrid = GetPropType<TypeTag, Properties::EquilGrid>;
using GridView = GetPropType<TypeTag, Properties::GridView>;
using TransmissibilityType = EclTransmissibility<Grid, GridView, ElementMapper, Scalar>;
static const int dimensionworld = Grid::dimensionworld;
private:
typedef Dune::CartesianIndexMapper<Grid> CartesianIndexMapper;
@@ -138,7 +139,6 @@ public:
this->eclState(), this->parallelWells_);
#endif
this->allocCartMapper();
this->updateGridView_();
this->updateCartesianToCompressedMapping_();
this->updateCellDepths_();
@@ -149,6 +149,18 @@ public:
#endif
}
/*!
* \brief Get function to query cell centroids for a distributed grid.
*
* Currently this only non-empty for a loadbalanced CpGrid.
* It is a function return the centroid for the given element
* index.
*/
std::function<std::array<double,dimensionworld>(int)>
cellCentroids() const
{
return this->cellCentroids_(this->cartesianIndexMapper());
}
protected:
void createGrids_()

14
ebos/eclgenericcpgridvanguard.cc
View File

@@ -92,7 +92,6 @@ void EclGenericCpGridVanguard<ElementMapper,GridView,Scalar>::doLoadBalance_(Dun
// its edge weights. since this is (kind of) a layering violation and
// transmissibilities are relatively expensive to compute, we only do it if
// more than a single process is involved in the simulation.
cartesianIndexMapper_.reset(new CartesianIndexMapper(*grid_));
if (grid_->size(0))
{
this->allocTrans();
@@ -175,8 +174,6 @@ void EclGenericCpGridVanguard<ElementMapper,GridView,Scalar>::doLoadBalance_(Dun
}
grid_->switchToDistributedView();
cartesianIndexMapper_.reset();
// Calling Schedule::filterConnections would remove any perforated
// cells that exist only on other ranks even in the case of distributed wells
// But we need all connections to figure out the first cell of a well (e.g. for
@@ -210,12 +207,6 @@ void EclGenericCpGridVanguard<ElementMapper,GridView,Scalar>::distributeFieldPro
}
#endif
template<class ElementMapper, class GridView, class Scalar>
void EclGenericCpGridVanguard<ElementMapper,GridView,Scalar>::allocCartMapper()
{
this->cartesianIndexMapper_.reset(new CartesianIndexMapper(this->grid()));
}
template<class ElementMapper, class GridView, class Scalar>
void EclGenericCpGridVanguard<ElementMapper,GridView,Scalar>::doCreateGrids_(EclipseState& eclState)
{
@@ -258,6 +249,9 @@ void EclGenericCpGridVanguard<ElementMapper,GridView,Scalar>::doCreateGrids_(Ecl
}
}
cartesianIndexMapper_ = std::make_unique<CartesianIndexMapper>(*grid_);
#if HAVE_MPI
{
const bool has_numerical_aquifer = eclState.aquifer().hasNumericalAquifer();
@@ -285,7 +279,7 @@ void EclGenericCpGridVanguard<ElementMapper,GridView,Scalar>::doCreateGrids_(Ecl
if (mpiRank == 0)
{
equilGrid_.reset(new Dune::CpGrid(*grid_));
equilCartesianIndexMapper_.reset(new CartesianIndexMapper(*equilGrid_));
equilCartesianIndexMapper_ = std::make_unique<CartesianIndexMapper>(*equilGrid_);
std::vector<int> actnum = UgGridHelpers::createACTNUM(*grid_);
auto &field_props = eclState.fieldProps();

17
ebos/eclpolyhedralgridvanguard.hh
View File

@@ -92,7 +92,7 @@ private:
typedef Grid* GridPointer;
typedef EquilGrid* EquilGridPointer;
typedef Dune::CartesianIndexMapper<Grid> CartesianIndexMapper;
typedef CartesianIndexMapper* CartesianIndexMapperPointer;
typedef std::unique_ptr<CartesianIndexMapper> CartesianIndexMapperPointer;
public:
EclPolyhedralGridVanguard(Simulator& simulator)
@@ -104,7 +104,6 @@ public:
~EclPolyhedralGridVanguard()
{
delete cartesianIndexMapper_;
delete grid_;
}
@@ -184,11 +183,23 @@ public:
return simulator_.problem().eclTransmissibilities();
}
/*!
* \brief Get function to query cell centroids for a distributed grid.
*
* Currently this only non-empty for a loadbalanced CpGrid.
* It is a function return the centroid for the given element
* index.
*/
std::function<std::array<double,dimensionworld>(int)>
cellCentroids() const
{
return this->cellCentroids_(this->cartesianIndexMapper());
}
protected:
void createGrids_()
{
grid_ = new Grid(this->eclState().getInputGrid(), this->eclState().fieldProps().porv(true));
cartesianIndexMapper_ = new CartesianIndexMapper(*grid_);
cartesianIndexMapper_ = std::make_unique<CartesianIndexMapper>(*grid_);
this->updateCartesianToCompressedMapping_();
this->updateCellDepths_();
}

12
ebos/ecltransmissibility.cc
View File

@@ -88,7 +88,7 @@ EclTransmissibility(const EclipseState& eclState,
const GridView& gridView,
const Dune::CartesianIndexMapper<Grid>& cartMapper,
const Grid& grid,
const std::vector<double>& centroids,
std::function<std::array<double,dimWorld>(int)> centroids,
bool enableEnergy,
bool enableDiffusivity)
: eclState_(eclState)
@@ -174,15 +174,7 @@ update(bool global)
// compute the axis specific "centroids" used for the transmissibilities. for
// consistency with the flow simulator, we use the element centers as
// computed by opm-parser's Opm::EclipseGrid class for all axes.
std::array<double, 3> centroid;
if (gridView_.comm().rank() == 0) {
const auto& eclGrid = eclState_.getInputGrid();
unsigned cartesianCellIdx = cartMapper_.cartesianIndex(elemIdx);
centroid = eclGrid.getCellCenter(cartesianCellIdx);
} else
std::copy(centroids_.begin() + centroidIdx * dimWorld,
centroids_.begin() + (centroidIdx + 1) * dimWorld,
centroid.begin());
std::array<double, dimWorld> centroid = centroids_(elemIdx);
for (unsigned axisIdx = 0; axisIdx < dimWorld; ++axisIdx)
for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)

5
ebos/ecltransmissibility.hh
View File

@@ -38,6 +38,7 @@
#include <tuple>
#include <vector>
#include <unordered_map>
#include <functional>
namespace Opm {
@@ -58,7 +59,7 @@ public:
const GridView& gridView,
const Dune::CartesianIndexMapper<Grid>& cartMapper,
const Grid& grid,
const std::vector<double>& centroids,
std::function<std::array<double,dimWorld>(int)> centroids,
bool enableEnergy,
bool enableDiffusivity);
@@ -236,7 +237,7 @@ protected:
const GridView& gridView_;
const Dune::CartesianIndexMapper<Grid>& cartMapper_;
const Grid& grid_;
const std::vector<double>& centroids_;
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_;