diff --git a/dumux/implicit/box/stencils.hh b/dumux/implicit/box/stencils.hh new file mode 100644 index 0000000000000000000000000000000000000000..19ef1b63e713ca750ce7323d0ea003dc21872e79 --- /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 e4fcc3e503815bbcbfa0202d9229fda439ecd0bc..ad6f15b3e37b75c54782b34c0b6ac93096159a39 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 ae263ee82236848fb218e280d1c09c97ad9a3216..edbd34c329acba9429a936ac3cea8088a8e9f9b7 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 978d76a50d48e5ce5a974648fa861e59304ef59c..41e22f039f8d15d6a7f5d14afcfe8c26109c37e7 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 c63c3b9da61a54303be462bf412cdd40c7724fbe..213858b2bf66fad23b6fb7455f268c283bba5a4f 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 83768e0d58cc80318effe1e1f5f7f7df2cedc3c6..0000000000000000000000000000000000000000 --- 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