Commit 43dbcd24 authored by Ned Coltman's avatar Ned Coltman
Browse files

[discretization][box] replace all ReferenceElements::general calls with referenceElement

parent a614e125
......@@ -26,7 +26,6 @@
#define DUMUX_DISCRETIZATION_BOX_GEOMETRY_HELPER_HH
#include <array>
#include <dune/geometry/referenceelements.hh>
#include <dumux/common/math.hh>
......@@ -128,7 +127,6 @@ class BoxGeometryHelper<GridView, 2, ScvType, ScvfType>
static constexpr auto dim = GridView::dimension;
static constexpr auto dimWorld = GridView::dimensionworld;
using ReferenceElements = typename Dune::ReferenceElements<Scalar, dim>;
//! the maximum number of helper points used to construct the geometries
//! Using a statically sized point array is much faster than dynamic allocation
......@@ -140,7 +138,7 @@ public:
: elementGeometry_(geometry)
, corners_(geometry.corners())
{
const auto referenceElement = ReferenceElements::general(geometry.type());
const auto refElement = referenceElement(geometry);
// the element center
p_[0] = geometry.center();
......@@ -150,8 +148,8 @@ public:
p_[i+1] = geometry.corner(i);
// face midpoints
for (int i = 0; i < referenceElement.size(1); ++i)
p_[i+corners_+1] = geometry.global(referenceElement.position(i, 1));
for (int i = 0; i < refElement.size(1); ++i)
p_[i+corners_+1] = geometry.global(refElement.position(i, 1));
}
//! Create a vector with the scv corners
......@@ -250,9 +248,9 @@ public:
const typename Intersection::Geometry& isGeom,
unsigned int indexInIntersection) const
{
const auto referenceElement = ReferenceElements::general(elementGeometry_.type());
const auto refElement = referenceElement(elementGeometry_);
const auto vIdxLocal = referenceElement.subEntity(is.indexInInside(), 1, indexInIntersection, dim);
const auto vIdxLocal = refElement.subEntity(is.indexInInside(), 1, indexInIntersection, dim);
if (indexInIntersection == 0)
return ScvfCornerStorage({p_[vIdxLocal+1], isGeom.center()});
else if (indexInIntersection == 1)
......@@ -348,8 +346,6 @@ class BoxGeometryHelper<GridView, 3, ScvType, ScvfType>
static constexpr auto dim = GridView::dimension;
static constexpr auto dimWorld = GridView::dimensionworld;
using ReferenceElements = typename Dune::ReferenceElements<Scalar, dim>;
using FaceReferenceElements = typename Dune::ReferenceElements<Scalar, dim-1>;
//! the maximum number of helper points used to construct the geometries
//! Using a statically sized point array is much faster than dynamic allocation
......@@ -359,7 +355,7 @@ public:
: elementGeometry_(geometry)
, corners_(geometry.corners())
{
const auto referenceElement = ReferenceElements::general(geometry.type());
const auto refElement = referenceElement(geometry);
// the element center
p_[0] = geometry.center();
......@@ -369,12 +365,12 @@ public:
p_[i+1] = geometry.corner(i);
// edge midpoints
for (int i = 0; i < referenceElement.size(dim-1); ++i)
p_[i+corners_+1] = geometry.global(referenceElement.position(i, dim-1));
for (int i = 0; i < refElement.size(dim-1); ++i)
p_[i+corners_+1] = geometry.global(refElement.position(i, dim-1));
// face midpoints
for (int i = 0; i < referenceElement.size(1); ++i)
p_[i+corners_+1+referenceElement.size(dim-1)] = geometry.global(referenceElement.position(i, 1));
for (int i = 0; i < refElement.size(1); ++i)
p_[i+corners_+1+refElement.size(dim-1)] = geometry.global(refElement.position(i, 1));
}
//! Create a vector with the scv corners
......@@ -552,8 +548,8 @@ public:
const typename Intersection::Geometry& geometry,
unsigned int indexInIntersection) const
{
const auto referenceElement = ReferenceElements::general(elementGeometry_.type());
const auto faceRefElem = FaceReferenceElements::general(geometry.type());
const auto refElement = referenceElement(elementGeometry_);
const auto faceRefElem = referenceElement(geometry);
GlobalPosition pi[9];
auto corners = geometry.corners();
......@@ -565,14 +561,14 @@ public:
const auto idxInInside = is.indexInInside();
for (int i = 0; i < corners; ++i)
{
const auto vIdxLocal = referenceElement.subEntity(idxInInside, 1, i, dim);
const auto vIdxLocal = refElement.subEntity(idxInInside, 1, i, dim);
pi[i+1] = elementGeometry_.corner(vIdxLocal);
}
// edge midpoints
for (int i = 0; i < faceRefElem.size(1); ++i)
{
const auto edgeIdxLocal = referenceElement.subEntity(idxInInside, 1, i, dim-1);
const auto edgeIdxLocal = refElement.subEntity(idxInInside, 1, i, dim-1);
pi[i+corners+1] = p_[edgeIdxLocal+corners_+1];
}
......
......@@ -27,7 +27,6 @@
#define DUMUX_DISCRETIZATION_BOX_FV_ELEMENT_GEOMETRY_HH
#include <dune/geometry/type.hh>
#include <dune/geometry/referenceelements.hh>
#include <dune/localfunctions/lagrange/pqkfactory.hh>
#include <dumux/common/indextraits.hh>
......@@ -58,7 +57,6 @@ class BoxFVElementGeometry<GG, true>
using LocalIndexType = typename IndexTraits<GridView>::LocalIndex;
using CoordScalar = typename GridView::ctype;
using FeLocalBasis = typename GG::FeCache::FiniteElementType::Traits::LocalBasisType;
using ReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim>;
public:
//! export the element type
using Element = typename GridView::template Codim<0>::Entity;
......@@ -174,7 +172,6 @@ class BoxFVElementGeometry<GG, false>
using LocalIndexType = typename IndexTraits<GridView>::LocalIndex;
using CoordScalar = typename GridView::ctype;
using FeLocalBasis = typename GG::FeCache::FiniteElementType::Traits::LocalBasisType;
using ReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim>;
using GeometryHelper = BoxGeometryHelper<GridView, dim,
typename GG::SubControlVolume,
......@@ -283,7 +280,7 @@ private:
// get the element geometry
auto elementGeometry = element.geometry();
elemGeometryType_ = elementGeometry.type();
const auto referenceElement = ReferenceElements::general(elemGeometryType_);
const auto refElement = referenceElement(elementGeometry);
// get the sub control volume geometries of this element
GeometryHelper geometryHelper(elementGeometry);
......@@ -310,8 +307,8 @@ private:
for (; scvfLocalIdx < numInnerScvf; ++scvfLocalIdx)
{
// find the local scv indices this scvf is connected to
std::vector<LocalIndexType> localScvIndices({static_cast<LocalIndexType>(referenceElement.subEntity(scvfLocalIdx, dim-1, 0, dim)),
static_cast<LocalIndexType>(referenceElement.subEntity(scvfLocalIdx, dim-1, 1, dim))});
std::vector<LocalIndexType> localScvIndices({static_cast<LocalIndexType>(refElement.subEntity(scvfLocalIdx, dim-1, 0, dim)),
static_cast<LocalIndexType>(refElement.subEntity(scvfLocalIdx, dim-1, 1, dim))});
scvfs_[scvfLocalIdx] = SubControlVolumeFace(geometryHelper,
element,
......@@ -332,7 +329,7 @@ private:
for (unsigned int isScvfLocalIdx = 0; isScvfLocalIdx < isGeometry.corners(); ++isScvfLocalIdx)
{
// find the scv this scvf is connected to
const LocalIndexType insideScvIdx = static_cast<LocalIndexType>(referenceElement.subEntity(intersection.indexInInside(), 1, isScvfLocalIdx, dim));
const LocalIndexType insideScvIdx = static_cast<LocalIndexType>(refElement.subEntity(intersection.indexInInside(), 1, isScvfLocalIdx, dim));
std::vector<LocalIndexType> localScvIndices = {insideScvIdx, insideScvIdx};
scvfs_.emplace_back(geometryHelper,
......
......@@ -28,7 +28,6 @@
#include <unordered_map>
#include <dune/geometry/referenceelements.hh>
#include <dune/localfunctions/lagrange/pqkfactory.hh>
#include <dumux/discretization/method.hh>
......@@ -91,8 +90,6 @@ class BoxFVGridGeometry<Scalar, GV, true, Traits>
static const int dim = GV::dimension;
static const int dimWorld = GV::dimensionworld;
using ReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim>;
using GeometryHelper = BoxGeometryHelper<GV, dim,
typename Traits::SubControlVolume,
typename Traits::SubControlVolumeFace>;
......@@ -171,7 +168,7 @@ public:
// get the element geometry
auto elementGeometry = element.geometry();
const auto referenceElement = ReferenceElements::general(elementGeometry.type());
const auto refElement = referenceElement(elementGeometry);
// instantiate the geometry helper
GeometryHelper geometryHelper(elementGeometry);
......@@ -194,8 +191,8 @@ public:
for (; scvfLocalIdx < element.subEntities(dim-1); ++scvfLocalIdx)
{
// find the global and local scv indices this scvf is belonging to
std::vector<LocalIndexType> localScvIndices({static_cast<LocalIndexType>(referenceElement.subEntity(scvfLocalIdx, dim-1, 0, dim)),
static_cast<LocalIndexType>(referenceElement.subEntity(scvfLocalIdx, dim-1, 1, dim))});
std::vector<LocalIndexType> localScvIndices({static_cast<LocalIndexType>(refElement.subEntity(scvfLocalIdx, dim-1, 0, dim)),
static_cast<LocalIndexType>(refElement.subEntity(scvfLocalIdx, dim-1, 1, dim))});
scvfs_[eIdx][scvfLocalIdx] = SubControlVolumeFace(geometryHelper,
element,
......@@ -220,7 +217,7 @@ public:
for (unsigned int isScvfLocalIdx = 0; isScvfLocalIdx < isGeometry.corners(); ++isScvfLocalIdx)
{
// find the scvs this scvf is belonging to
const LocalIndexType insideScvIdx = static_cast<LocalIndexType>(referenceElement.subEntity(intersection.indexInInside(), 1, isScvfLocalIdx, dim));
const LocalIndexType insideScvIdx = static_cast<LocalIndexType>(refElement.subEntity(intersection.indexInInside(), 1, isScvfLocalIdx, dim));
std::vector<LocalIndexType> localScvIndices = {insideScvIdx, insideScvIdx};
scvfs_[eIdx].emplace_back(geometryHelper,
......@@ -238,10 +235,10 @@ public:
// add all vertices on the intersection to the set of
// boundary vertices
const auto fIdx = intersection.indexInInside();
const auto numFaceVerts = referenceElement.size(fIdx, 1, dim);
const auto numFaceVerts = refElement.size(fIdx, 1, dim);
for (int localVIdx = 0; localVIdx < numFaceVerts; ++localVIdx)
{
const auto vIdx = referenceElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdx = refElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdxGlobal = this->vertexMapper().subIndex(element, vIdx, dim);
boundaryDofIndices_[vIdxGlobal] = true;
}
......@@ -254,11 +251,11 @@ public:
// find the mapped periodic vertex of all vertices on periodic boundaries
const auto fIdx = intersection.indexInInside();
const auto numFaceVerts = referenceElement.size(fIdx, 1, dim);
const auto numFaceVerts = refElement.size(fIdx, 1, dim);
const auto eps = 1e-7*(elementGeometry.corner(1) - elementGeometry.corner(0)).two_norm();
for (int localVIdx = 0; localVIdx < numFaceVerts; ++localVIdx)
{
const auto vIdx = referenceElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdx = refElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdxGlobal = this->vertexMapper().subIndex(element, vIdx, dim);
const auto vPos = elementGeometry.corner(vIdx);
......@@ -270,10 +267,10 @@ public:
if (isOutside.boundary() && isOutside.neighbor())
{
const auto fIdxOutside = isOutside.indexInInside();
const auto numFaceVertsOutside = referenceElement.size(fIdxOutside, 1, dim);
const auto numFaceVertsOutside = refElement.size(fIdxOutside, 1, dim);
for (int localVIdxOutside = 0; localVIdxOutside < numFaceVertsOutside; ++localVIdxOutside)
{
const auto vIdxOutside = referenceElement.subEntity(fIdxOutside, 1, localVIdxOutside, dim);
const auto vIdxOutside = refElement.subEntity(fIdxOutside, 1, localVIdxOutside, dim);
const auto vPosOutside = outsideGeometry.corner(vIdxOutside);
const auto shift = std::abs((this->bBoxMax()-this->bBoxMin())*intersection.centerUnitOuterNormal());
if (std::abs((vPosOutside-vPos).two_norm() - shift) < eps)
......@@ -363,8 +360,6 @@ class BoxFVGridGeometry<Scalar, GV, false, Traits>
using Element = typename GV::template Codim<0>::Entity;
using CoordScalar = typename GV::ctype;
using ReferenceElements = typename Dune::ReferenceElements<CoordScalar, dim>;
public:
//! export discretization method
static constexpr DiscretizationMethod discMethod = DiscretizationMethod::box;
......@@ -427,7 +422,7 @@ public:
numScvf_ += element.subEntities(dim-1);
const auto elementGeometry = element.geometry();
const auto referenceElement = ReferenceElements::general(elementGeometry.type());
const auto refElement = referenceElement(elementGeometry);
// store the sub control volume face indices on the domain boundary
for (const auto& intersection : intersections(this->gridView(), element))
......@@ -441,10 +436,10 @@ public:
// add all vertices on the intersection to the set of
// boundary vertices
const auto fIdx = intersection.indexInInside();
const auto numFaceVerts = referenceElement.size(fIdx, 1, dim);
const auto numFaceVerts = refElement.size(fIdx, 1, dim);
for (int localVIdx = 0; localVIdx < numFaceVerts; ++localVIdx)
{
const auto vIdx = referenceElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdx = refElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdxGlobal = this->vertexMapper().subIndex(element, vIdx, dim);
boundaryDofIndices_[vIdxGlobal] = true;
}
......@@ -457,11 +452,11 @@ public:
// find the mapped periodic vertex of all vertices on periodic boundaries
const auto fIdx = intersection.indexInInside();
const auto numFaceVerts = referenceElement.size(fIdx, 1, dim);
const auto numFaceVerts = refElement.size(fIdx, 1, dim);
const auto eps = 1e-7*(elementGeometry.corner(1) - elementGeometry.corner(0)).two_norm();
for (int localVIdx = 0; localVIdx < numFaceVerts; ++localVIdx)
{
const auto vIdx = referenceElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdx = refElement.subEntity(fIdx, 1, localVIdx, dim);
const auto vIdxGlobal = this->vertexMapper().subIndex(element, vIdx, dim);
const auto vPos = elementGeometry.corner(vIdx);
......@@ -473,10 +468,10 @@ public:
if (isOutside.boundary() && isOutside.neighbor())
{
const auto fIdxOutside = isOutside.indexInInside();
const auto numFaceVertsOutside = referenceElement.size(fIdxOutside, 1, dim);
const auto numFaceVertsOutside = refElement.size(fIdxOutside, 1, dim);
for (int localVIdxOutside = 0; localVIdxOutside < numFaceVertsOutside; ++localVIdxOutside)
{
const auto vIdxOutside = referenceElement.subEntity(fIdxOutside, 1, localVIdxOutside, dim);
const auto vIdxOutside = refElement.subEntity(fIdxOutside, 1, localVIdxOutside, dim);
const auto vPosOutside = outsideGeometry.corner(vIdxOutside);
const auto shift = std::abs((this->bBoxMax()-this->bBoxMin())*intersection.centerUnitOuterNormal());
if (std::abs((vPosOutside-vPos).two_norm() - shift) < eps)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment