From 53751a567598f21d401b1e56a02714e48bc6c6ac Mon Sep 17 00:00:00 2001
From: Timo Koch <timo.koch@iws.uni-stuttgart.de>
Date: Mon, 21 Mar 2016 01:56:42 +0100
Subject: [PATCH] [stencils] Generalize stencils
There are now stencils related to elements and stencils related to vertices.
The element stencil for box e.g. contains all dof related to an element
whereas the vertex stencil contains all dofs related to a vertex (==dof).
---
dumux/implicit/box/stencils.hh | 159 ++++++++++++++++++++++++
dumux/implicit/cellcentered/stencils.hh | 44 ++++++-
dumux/implicit/model.hh | 15 +--
dumux/implicit/properties.hh | 1 -
dumux/implicit/propertydefaults.hh | 4 -
dumux/implicit/stencilsvector.hh | 54 --------
6 files changed, 206 insertions(+), 71 deletions(-)
create mode 100644 dumux/implicit/box/stencils.hh
delete mode 100644 dumux/implicit/stencilsvector.hh
diff --git a/dumux/implicit/box/stencils.hh b/dumux/implicit/box/stencils.hh
new file mode 100644
index 0000000000..19ef1b63e7
--- /dev/null
+++ b/dumux/implicit/box/stencils.hh
@@ -0,0 +1,159 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ * See the file COPYING for full copying permissions. *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 2 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ *****************************************************************************/
+/*!
+ * \file
+ * \brief Implements the notion of stencils for vertex-centered models
+ */
+#ifndef DUMUX_BOX_STENCILS_HH
+#define DUMUX_BOX_STENCILS_HH
+
+#include <set>
+#include <dumux/implicit/box/properties.hh>
+
+namespace Dumux
+{
+//forward declaration
+template<class TypeTag>
+class BoxStencilsVector;
+
+/*!
+ * \brief Element-related stencils
+ */
+template<class TypeTag>
+class BoxElementStencils
+{
+ using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+ using IndexType = typename GridView::IndexSet::IndexType;
+ using Element = typename GridView::template Codim<0>::Entity;
+ // TODO a separate stencil class all stencils can derive from?
+ using Stencil = std::set<IndexType>;
+
+ static const int dim = GridView::dimension;
+public:
+ void update(const Problem& problem, const Element& element)
+ {
+ for(int vIdxLocal = 0; vIdxLocal < element.subEntities(dim); ++vIdxLocal)
+ elementStencil_.insert(problem.vertexMapper().subIndex(element, vIdxLocal, dim));
+ }
+
+ //! The full element stencil (all element this element is interacting with)
+ const Stencil& elementStencil() const
+ {
+ return elementStencil_;
+ }
+
+private:
+ Stencil elementStencil_;
+};
+
+/*!
+ * \brief Vertex-related stencils
+ */
+template<class TypeTag>
+class BoxVertexStencils
+{
+ friend class BoxStencilsVector<TypeTag>;
+ using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+ using IndexType = typename GridView::IndexSet::IndexType;
+ using Element = typename GridView::template Codim<0>::Entity;
+ // TODO a separate stencil class all stencils can derive from?
+ using Stencil = std::set<IndexType>;
+public:
+ //! The full vertex stencil (all vertices this vertex is interacting with)
+ const Stencil& vertexStencil() const
+ {
+ return vertexStencil_;
+ }
+
+private:
+ //! The full vertex stencil (all vertices this vertex is interacting with)
+ Stencil& vertexStencil()
+ {
+ return vertexStencil_;
+ }
+
+ Stencil vertexStencil_;
+};
+
+/*!
+ * \ingroup BoxModel
+ * \brief The global stencil container class
+ */
+template<class TypeTag>
+class BoxStencilsVector
+{
+ using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+
+ static const int dim = GridView::dimension;
+
+public:
+ void update(const Problem& problem)
+ {
+ problemPtr_ = &problem;
+ elementStencils_.resize(problem.gridView().size(0));
+ vertexStencils_.resize(problem.gridView().size(dim));
+ for (const auto& element : elements(problem.gridView()))
+ {
+ auto eIdx = problem.elementMapper().index(element);
+ elementStencils_[eIdx].update(problem, element);
+
+ for (int vIdxLocalI = 0; vIdxLocalI < element.subEntities(dim); ++vIdxLocalI)
+ {
+ auto globalI = problem.vertexMapper().subIndex(element, vIdxLocalI, dim);
+ for (int vIdxLocalJ = vIdxLocalI; vIdxLocalJ < element.subEntities(dim); ++vIdxLocalJ)
+ {
+ auto globalJ = problem.vertexMapper().subIndex(element, vIdxLocalJ, dim);
+
+ // make sure that vertex j is in the neighbor set
+ // of vertex i and vice-versa
+ vertexStencils_[globalI].vertexStencil().insert(globalJ);
+ vertexStencils_[globalJ].vertexStencil().insert(globalI);
+ }
+ }
+ }
+ }
+
+ //! overload for elements
+ template <class Entity>
+ typename std::enable_if<Entity::codimension == 0, const BoxElementStencils<TypeTag>&>::type
+ get(const Entity& entity) const
+ {
+ return elementStencils_[problemPtr_->elementMapper().index(entity)];
+ }
+
+ //! overload for vertices
+ template <class Entity>
+ typename std::enable_if<Entity::codimension == Entity::dimension, const BoxVertexStencils<TypeTag>&>::type
+ get(const Entity& entity) const
+ {
+ return vertexStencils_[problemPtr_->vertexMapper().index(entity)];
+ }
+
+private:
+ std::vector<BoxElementStencils<TypeTag>> elementStencils_;
+ std::vector<BoxVertexStencils<TypeTag>> vertexStencils_;
+ const Problem* problemPtr_;
+};
+
+} // end namespace
+
+#endif
diff --git a/dumux/implicit/cellcentered/stencils.hh b/dumux/implicit/cellcentered/stencils.hh
index e4fcc3e503..ad6f15b3e3 100644
--- a/dumux/implicit/cellcentered/stencils.hh
+++ b/dumux/implicit/cellcentered/stencils.hh
@@ -18,7 +18,7 @@
*****************************************************************************/
/*!
* \file
- * \brief Implements the notion of stencils
+ * \brief Implements the notion of stencils for cell-centered models
*/
#ifndef DUMUX_CC_STENCILS_HH
#define DUMUX_CC_STENCILS_HH
@@ -28,13 +28,12 @@
namespace Dumux
{
+
/*!
- * \ingroup ImplicitModel
- * \brief Base class for a sub control volume, i.e a part of the control
- * volume we are making the balance for.
+ * \brief Element-related stencils
*/
template<class TypeTag>
-class CCStencils
+class CCElementStencils
{
using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
@@ -71,6 +70,41 @@ private:
Stencil neighborStencil_;
};
+/*!
+ * \ingroup CCModel
+ * \brief The global stencil container class
+ */
+template<class TypeTag>
+class CCStencilsVector
+{
+ using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+
+public:
+ void update(const Problem& problem)
+ {
+ problemPtr_ = &problem;
+ elementStencils_.resize(problem.gridView().size(0));
+ for (const auto& element : elements(problem.gridView()))
+ {
+ auto eIdx = problem.elementMapper().index(element);
+ elementStencils_[eIdx].update(problem, element);
+ }
+ }
+
+ //! overload for elements
+ template <class Entity>
+ typename std::enable_if<Entity::codimension == 0, const CCElementStencils<TypeTag>&>::type
+ get(const Entity& entity) const
+ {
+ return elementStencils_[problemPtr_->elementMapper().index(entity)];
+ }
+
+private:
+ std::vector<CCElementStencils<TypeTag>> elementStencils_;
+ const Problem* problemPtr_;
+};
+
} // end namespace
#endif
diff --git a/dumux/implicit/model.hh b/dumux/implicit/model.hh
index ae263ee822..edbd34c329 100644
--- a/dumux/implicit/model.hh
+++ b/dumux/implicit/model.hh
@@ -59,7 +59,6 @@ class ImplicitModel
typedef typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace) SubControlVolumeFace;
typedef typename GET_PROP_TYPE(TypeTag, FluxVariables) FluxVariables;
typedef typename GET_PROP_TYPE(TypeTag, FluxVariablesVector) FluxVariablesVector;
- typedef typename GET_PROP_TYPE(TypeTag, Stencils) Stencils;
typedef typename GET_PROP_TYPE(TypeTag, StencilsVector) StencilsVector;
enum {
@@ -778,12 +777,6 @@ public:
const FVElementGeometryVector& fvGeometries() const
{ return *fvGeometries_; }
- const Stencils& stencils(const Element& element) const
- { return stencilsVector_[elementMapper().index(element)]; }
-
- const Stencils& stencils(unsigned int eIdx) const
- { return stencilsVector_[eIdx]; }
-
const FVElementGeometry& fvGeometries(const Element& element) const
{ return fvGeometries_->fvGeometry(elementMapper().index(element)); }
@@ -976,6 +969,14 @@ protected:
Dune::BlockVector<Dune::FieldVector<Scalar, 1> > boxVolume_;
+public:
+ //! Get stencils related to an entity of specified codimension
+ // Cell-centered discretizations typically only implement element stencils
+ // Vertex-centered discretizations both vertex and element stencils
+ template <class Entity>
+ auto stencils(const Entity& entity) const -> decltype(stencilsVector_.get(entity))
+ { return stencilsVector_.get(entity); }
+
private:
/*!
* \brief Returns whether messages should be printed
diff --git a/dumux/implicit/properties.hh b/dumux/implicit/properties.hh
index 978d76a50d..41e22f039f 100644
--- a/dumux/implicit/properties.hh
+++ b/dumux/implicit/properties.hh
@@ -87,7 +87,6 @@ NEW_PROP_TAG(EnableEnergyBalance); //! Specifies if the model solves an energy e
NEW_PROP_TAG(HeatConductionType); //! The type for the calculation of the heat conduction fluxes
// stencils
-NEW_PROP_TAG(Stencils); //! The stencils container type
NEW_PROP_TAG(StencilsVector); //! The type of the global vector of stencils per element
// high level simulation control
diff --git a/dumux/implicit/propertydefaults.hh b/dumux/implicit/propertydefaults.hh
index c63c3b9da6..213858b2bf 100644
--- a/dumux/implicit/propertydefaults.hh
+++ b/dumux/implicit/propertydefaults.hh
@@ -44,7 +44,6 @@
#include "fluxvariables.hh"
#include "fluxvariablesvector.hh"
#include "fvelementgeometry.hh"
-#include "stencilsvector.hh"
namespace Dumux {
@@ -99,9 +98,6 @@ SET_TYPE_PROP(ImplicitBase, VolumeVariables, ImplicitVolumeVariables<TypeTag>);
//! The global volume variables vector class
SET_TYPE_PROP(ImplicitBase, VolumeVariablesVector, Dumux::VolumeVariablesVector<TypeTag>);
-//! The global volume variables vector class
-SET_TYPE_PROP(ImplicitBase, StencilsVector, Dumux::StencilsVector<TypeTag>);
-
//! The global volume variables vector class
SET_TYPE_PROP(ImplicitBase, FluxVariablesVector, Dumux::FluxVariablesVector<TypeTag>);
diff --git a/dumux/implicit/stencilsvector.hh b/dumux/implicit/stencilsvector.hh
deleted file mode 100644
index 83768e0d58..0000000000
--- a/dumux/implicit/stencilsvector.hh
+++ /dev/null
@@ -1,54 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * See the file COPYING for full copying permissions. *
- * *
- * This program is free software: you can redistribute it and/or modify *
- * it under the terms of the GNU General Public License as published by *
- * the Free Software Foundation, either version 2 of the License, or *
- * (at your option) any later version. *
- * *
- * This program is distributed in the hope that it will be useful, *
- * but WITHOUT ANY WARRANTY; without even the implied warranty of *
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
- * GNU General Public License for more details. *
- * *
- * You should have received a copy of the GNU General Public License *
- * along with this program. If not, see <http://www.gnu.org/licenses/>. *
- *****************************************************************************/
-/*!
- * \file
- * \brief Implements the notion of stencils
- */
-#ifndef DUMUX_STENCILSVECTOR_HH
-#define DUMUX_STENCILSVECTOR_HH
-
-#include <dumux/implicit/properties.hh>
-
-namespace Dumux
-{
-/*!
- * \ingroup ImplicitModel
- * \brief Base class for a sub control volume, i.e a part of the control
- * volume we are making the balance for.
- */
-template<class TypeTag>
-class StencilsVector : public std::vector<typename GET_PROP_TYPE(TypeTag, Stencils)>
-{
- using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
- using Stencils = typename GET_PROP_TYPE(TypeTag, Stencils);
-public:
- void update(const Problem& problem)
- {
- this->resize(problem.gridView().size(0));
- for (const auto& element : elements(problem.gridView()))
- {
- auto eIdx = problem.elementMapper().index(element);
- (*this)[eIdx].update(problem, element);
- }
- }
-};
-
-} // end namespace
-
-#endif
--
GitLab