ebos: precalculate the depth of the center of an element

turns out that this is quite slow in Dune::CpGrid.
This commit is contained in:
Andreas Lauser 2016-10-25 15:54:14 +02:00 committed by Robert Kloefkorn
parent 171ccb0926
commit 8784c81a4e
2 changed files with 40 additions and 3 deletions

View File

@ -210,7 +210,6 @@ protected:
Valgrind::SetUndefined(*this);
const auto& problem = elemCtx.problem();
const auto& grid = elemCtx.simulator().gridManager().grid();
const auto& stencil = elemCtx.stencil(timeIdx);
const auto& scvf = stencil.interiorFace(scvfIdx);
@ -237,8 +236,8 @@ protected:
// grids). The "good" solution would be to take the Z coordinate of the element
// centers, but since ECL seems to like to be inconsistent on that front, it
// needs to be done like here...
Scalar zIn = grid.cellCenterDepth(I);
Scalar zEx = grid.cellCenterDepth(J);
Scalar zIn = problem.dofCenterDepth(elemCtx, interiorDofIdx_, timeIdx);
Scalar zEx = problem.dofCenterDepth(elemCtx, exteriorDofIdx_, timeIdx);
// the distances from the DOF's depths. (i.e., the additional depth of the
// exterior DOF)

View File

@ -319,6 +319,7 @@ public:
this->gravity_[dim - 1] = 9.80665;
initFluidSystem_();
updateElementDepths_();
readRockParameters_();
readMaterialParameters_();
transmissibilities_.finishInit();
@ -613,6 +614,19 @@ public:
return porosity_[globalSpaceIdx];
}
/*!
* \brief Returns the depth of an degree of freedom [m]
*
* For ECL problems this is defined as the average of the depth of an element and is
* thus slightly different from the depth of an element's centroid.
*/
template <class Context>
Scalar dofCenterDepth(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
{
unsigned globalSpaceIdx = context.globalSpaceIndex(spaceIdx, timeIdx);
return elementCenterDepth_[globalSpaceIdx];
}
/*!
* \copydoc BlackoilProblem::rockCompressibility
*/
@ -791,6 +805,29 @@ private:
static bool enableEclOutput_()
{ return EWOMS_GET_PARAM(TypeTag, bool, EnableEclOutput); }
void updateElementDepths_()
{
const auto& gridManager = this->simulator().gridManager();
const auto& grid = gridManager.grid();
const auto& gridView = gridManager.gridView();
const auto& elemMapper = this->elementMapper();;
int numElements = gridView.size(/*codim=*/0);
elementCenterDepth_.resize(numElements);
auto elemIt = gridView.template begin</*codim=*/0>();
const auto& elemEndIt = gridView.template end</*codim=*/0>();
for (; elemIt != elemEndIt; ++elemIt) {
#if DUNE_VERSION_NEWER(DUNE_COMMON, 2,4)
unsigned elemIdx = elemMapper.index(*elemIt);
#else
unsigned elemIdx = elemMapper.map(*elemIt);
#endif
elementCenterDepth_[elemIdx] = grid.cellCenterDepth(elemIdx);
}
}
void readRockParameters_()
{
auto deck = this->simulator().gridManager().deck();
@ -1186,6 +1223,7 @@ private:
}
std::vector<Scalar> porosity_;
std::vector<Scalar> elementCenterDepth_;
std::vector<DimMatrix> intrinsicPermeability_;
EclTransmissibility<TypeTag> transmissibilities_;