diff --git a/dune/functions/functionspacebases/CMakeLists.txt b/dune/functions/functionspacebases/CMakeLists.txt index 4f51de0711d3d5fe12b5f8ca46a79b4b4fd31aed..302b81b55e346dff5d7157d46e2727a76001120f 100644 --- a/dune/functions/functionspacebases/CMakeLists.txt +++ b/dune/functions/functionspacebases/CMakeLists.txt @@ -11,7 +11,7 @@ install(FILES defaultlocalview.hh defaultnodetorangemap.hh flatmultiindex.hh - flatvectorbackend.hh + flatvectorview.hh hierarchicvectorwrapper.hh interpolate.hh lagrangebasis.hh diff --git a/dune/functions/functionspacebases/flatvectorbackend.hh b/dune/functions/functionspacebases/flatvectorbackend.hh deleted file mode 100644 index 8dabee743d49349675ae6604652503de80228459..0000000000000000000000000000000000000000 --- a/dune/functions/functionspacebases/flatvectorbackend.hh +++ /dev/null @@ -1,82 +0,0 @@ -// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- -// vi: set et ts=4 sw=2 sts=2: -#ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORBACKEND_HH -#define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORBACKEND_HH - - -#include <dune/common/concept.hh> - -#include <dune/functions/functionspacebases/concepts.hh> - - - - -namespace Dune { -namespace Functions { - - - -template<class V> -struct FlatVectorBackend -{ - - template<class VV, class Index, - typename std::enable_if< models<Concept::HasIndexAccess, VV, Index>(), int>::type = 0> - static auto getEntry(VV&& v, const Index& i) - ->decltype(v[i]) - { - return v[i]; - } - - template<class VV, class Index, - typename std::enable_if< not models<Concept::HasIndexAccess, VV, Index>(), int>::type = 0> - static auto getEntry(VV&& v, const Index& i) - ->decltype(v) - { - return std::forward<VV>(v); - } - - template<class VV, - typename std::enable_if< models<Concept::HasSizeMethod, VV>(), int>::type = 0> - static auto size(VV&& v) - ->decltype(v.size()) - { - return v.size(); - } - - template<class VV, - typename std::enable_if< not models<Concept::HasSizeMethod, VV>(), int>::type = 0> - static std::size_t size(VV&& v) - { - return 1; - } - -}; - - - - - -template<class K, int n, int m> -struct FlatVectorBackend<typename Dune::FieldMatrix<K, n, m> > -{ - - template<class VV, class Index> - static auto getEntry(VV&& v, const Index& i) -> decltype(v[i/m][i%m]) - { - return v[i/m][i%m]; - } - - template<class VV> - static int size(VV&& v) - { - return n*m; - } -}; - - -} // namespace Dune::Functions -} // namespace Dune - - -#endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORBACKEND_HH diff --git a/dune/functions/functionspacebases/flatvectorview.hh b/dune/functions/functionspacebases/flatvectorview.hh new file mode 100644 index 0000000000000000000000000000000000000000..8a1afc015635c792decca4503dd0755480cf4f9b --- /dev/null +++ b/dune/functions/functionspacebases/flatvectorview.hh @@ -0,0 +1,205 @@ +// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- +// vi: set et ts=4 sw=2 sts=2: +#ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORVIEW_HH +#define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORVIEW_HH + + +#include <dune/common/concept.hh> + +#include <dune/functions/functionspacebases/concepts.hh> + + + + +namespace Dune { +namespace Functions { +namespace Impl { + + + +template<class V> +struct FlatVectorBackend +{ + + template<class VV, class Index, + typename std::enable_if< models<Concept::HasIndexAccess, VV, Index>(), int>::type = 0> + static auto getEntry(VV&& v, const Index& i) + ->decltype(v[i]) + { + return v[i]; + } + + template<class VV, class Index, + typename std::enable_if< not models<Concept::HasIndexAccess, VV, Index>(), int>::type = 0> + static auto getEntry(VV&& v, const Index& i) + ->decltype(v) + { + return std::forward<VV>(v); + } + + template<class VV, + typename std::enable_if< models<Concept::HasSizeMethod, VV>(), int>::type = 0> + static auto size(VV&& v) + ->decltype(v.size()) + { + return v.size(); + } + + template<class VV, + typename std::enable_if< not models<Concept::HasSizeMethod, VV>(), int>::type = 0> + static std::size_t size(VV&& v) + { + return 1; + } + +}; + + + + + +template<class K, int n, int m> +struct FlatVectorBackend<typename Dune::FieldMatrix<K, n, m> > +{ + + template<class VV, class Index> + static auto getEntry(VV&& v, const Index& i) -> decltype(v[i/m][i%m]) + { + return v[i/m][i%m]; + } + + template<class VV> + static int size(VV&& v) + { + return n*m; + } +}; + + + + +template<class T> +class FlatVectorView +{ + using Backend = FlatVectorBackend<std::decay_t<T>>; +public: + FlatVectorView(T& t) : + t_(&t) + {} + + auto size() const + { + return Backend::size(*t_); + } + + template<class Index> + decltype(auto) operator[](const Index& i) const + { + return Backend::getEntry(*t_, i); + } + + template<class Index> + decltype(auto) operator[](const Index& i) + { + return Backend::getEntry(*t_, i); + } + +private: + T* t_; +}; + + +template<class T> +class FlatVectorView<T&&> +{ + using Backend = FlatVectorBackend<std::decay_t<T>>; +public: + FlatVectorView(T&& t) : + t_(std::move(t)) + {} + + auto size() const + { + return Backend::size(t_); + } + + template<class Index> + decltype(auto) operator[](const Index& i) const + { + return Backend::getEntry(t_, i); + } + + template<class Index> + decltype(auto) operator[](const Index& i) + { + return Backend::getEntry(t_, i); + } + +private: + T t_; +}; + +} // namespace Impl + + + +/** + * \brief Create flat vector view of passed mutable container + * + * When passed a nested container, the resulting value is + * a flat-vector-like view object. It provides an operator[] + * method to access all entries of the underlying nested + * container using flat consecutive indices and a size() + * method to compute the corresponding total size. + * + * This method will create a view object storing a pointer + * to the passed mutable container. + */ +template<class T> +auto flatVectorView(T& t) +{ + return Impl::FlatVectorView<T>(t); +} + +/** + * \brief Create flat vector view of passed const container + * + * When passed a nested container, the resulting value is + * a flat-vector-like view object. It provides an operator[] + * method to access all entries of the underlying nested + * container using flat consecutive indices and a size() + * method to compute the corresponding total size. + * + * This method will create a view object storing a pointer + * to the passed const container. + */ +template<class T> +auto flatVectorView(const T& t) +{ + return Impl::FlatVectorView<const T>(t); +} + +/** + * \brief Create flat vector view of passed container temporary + * + * When passed a nested container, the resulting value is + * a flat-vector-like view object. It provides an operator[] + * method to access all entries of the underlying nested + * container using flat consecutive indices and a size() + * method to compute the corresponding total size. + * + * This method will create a 'view' object storing the + * provided temporary container by value. + */ +template<class T> +auto flatVectorView(T&& t) +{ + return Impl::FlatVectorView<T&&>(t); +} + + +} // namespace Dune::Functions +} // namespace Dune + + +#endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORVIEW_HH diff --git a/dune/functions/functionspacebases/interpolate.hh b/dune/functions/functionspacebases/interpolate.hh index f8f71b97ed3a415967d466e14656bf3ec977d28f..aa879981439c229438b1c868de784f867ee823d4 100644 --- a/dune/functions/functionspacebases/interpolate.hh +++ b/dune/functions/functionspacebases/interpolate.hh @@ -19,7 +19,7 @@ #include <dune/functions/functionspacebases/hierarchicvectorwrapper.hh> #include <dune/functions/functionspacebases/sizeinfo.hh> -#include <dune/functions/functionspacebases/flatvectorbackend.hh> +#include <dune/functions/functionspacebases/flatvectorview.hh> #include <dune/functions/functionspacebases/defaultnodetorangemap.hh> #include <dune/typetree/traversal.hh> @@ -84,9 +84,6 @@ public: using GlobalDomain = typename Element::Geometry::GlobalCoordinate; - using CoefficientBlock = typename std::decay<decltype(std::declval<HierarchicVector>()[std::declval<MultiIndex>()])>::type; - using BitVectorBlock = typename std::decay<decltype(std::declval<HierarchicBitVector>()[std::declval<MultiIndex>()])>::type; - LocalInterpolateVisitor(const B& basis, HV& coeff, const HBV& bitVector, const LF& localF, const LocalIndexSet& localIndexSet, const NodeToRangeEntry& nodeToRangeEntry) : vector_(coeff), localF_(localF), @@ -120,9 +117,7 @@ public: std::size_t j=0; auto localFj = [&](const LocalDomain& x){ const auto& y = localF_(x); - const auto& y_node = nodeToRangeEntry_(node, y); - using FunctionRange = typename std::decay<decltype(y_node)>::type; - return FlatVectorBackend<FunctionRange>::getEntry(y_node, j); + return flatVectorView(nodeToRangeEntry_(node, y))[j]; }; using FunctionFromCallable = typename Dune::Functions::FunctionFromCallable<FiniteElementRange(LocalDomain), decltype(localFj), FunctionBaseClass>; @@ -134,7 +129,7 @@ public: // We loop over j defined above and thus over the components of the // range type of localF_. - auto blockSize = FlatVectorBackend<CoefficientBlock>::size(vector_[localIndexSet_.index(0)]); + auto blockSize = flatVectorView(vector_[localIndexSet_.index(0)]).size(); for(j=0; j<blockSize; ++j) { @@ -144,13 +139,11 @@ public: for (size_t i=0; i<fe.localBasis().size(); ++i) { auto multiIndex = localIndexSet_.index(node.localIndex(i)); - const auto& bitVectorBlock = bitVector_[multiIndex]; - const auto& interpolateHere = FlatVectorBackend<BitVectorBlock>::getEntry(bitVectorBlock,j); - - if (interpolateHere) + auto bitVectorBlock = flatVectorView(bitVector_[multiIndex]); + if (bitVectorBlock[j]) { - auto&& vectorBlock = vector_[multiIndex]; - FlatVectorBackend<CoefficientBlock>::getEntry(vectorBlock, j) = interpolationCoefficients[i]; + auto vectorBlock = flatVectorView(vector_[multiIndex]); + vectorBlock[j] = interpolationCoefficients[i]; } } } @@ -180,7 +173,7 @@ protected: * * Notice that this will only work if the range type of f and * the block type of coeff are compatible and supported by - * FlatVectorBackend. + * flatVectorView. * * \param basis Global function space basis of discrete function space * \param treePath Tree path specifying the part of the ansatz tree to use @@ -262,7 +255,7 @@ void interpolateTree(const B& basis, const TypeTree::HybridTreePath<TreeIndices. * * Notice that this will only work if the range type of f and * the block type of coeff are compatible and supported by - * FlatVectorBackend. + * flatVectorView. * * \param basis Global function space basis of discrete function space * \param treePath Tree path specifying the part of the ansatz tree to use @@ -301,7 +294,7 @@ namespace Imp { * * Notice that this will only work if the range type of f and * the block type of coeff are compatible and supported by - * FlatVectorBackend. + * flatVectorView. * * \param basis Global function space basis of discrete function space * \param coeff Coefficient vector to represent the interpolation @@ -323,7 +316,7 @@ void interpolate(const B& basis, C&& coeff, const F& f, const BV& bitVector) * * Notice that this will only work if the range type of f and * the block type of coeff are compatible and supported by - * FlatVectorBackend. + * flatVectorView. * * This function will only work, if the local ansatz tree of * the basis is trivial, i.e., a single leaf node. @@ -346,7 +339,7 @@ void interpolate(const B& basis, C&& coeff, const F& f) * * Notice that this will only work if the range type of f and * the block type of corresponding coeff entries are compatible - * and supported by FlatVectorBackend. + * and supported by flatVectorView. * * \param basis Global function space basis of discrete function space * \param treePath Tree path specifying the part of the ansatz tree to use diff --git a/dune/functions/gridfunctions/discreteglobalbasisfunction.hh b/dune/functions/gridfunctions/discreteglobalbasisfunction.hh index 8e926d3e0fc24a644aab88bbbeffd0de6852db1d..807b02c0f59222e70df190e7077ae70476a951af 100644 --- a/dune/functions/gridfunctions/discreteglobalbasisfunction.hh +++ b/dune/functions/gridfunctions/discreteglobalbasisfunction.hh @@ -8,7 +8,7 @@ #include <dune/common/shared_ptr.hh> #include <dune/functions/functionspacebases/defaultnodetorangemap.hh> -#include <dune/functions/functionspacebases/flatvectorbackend.hh> +#include <dune/functions/functionspacebases/flatvectorview.hh> #include <dune/functions/gridfunctions/gridviewentityset.hh> #include <dune/functions/gridfunctions/gridfunction.hh> #include <dune/functions/common/treedata.hh> @@ -40,13 +40,13 @@ namespace Functions { * V the value of this basis function at the evaluation point. Notice * that both may be scalar, vector, matrix, or general container valued. * - * 2.Each entry of C is associated with a flat index j via FlatVectorBackend. + * 2.Each entry of C is associated with a flat index j via flatVectorView. * This is normally a lexicographic index. The total scalar dimension according * to those flat indices is dim(C). - * 3.Each entry of V is associated with a flat index k via FlatVectorBackend. + * 3.Each entry of V is associated with a flat index k via flatVectorView. * This is normally a lexicographic index. The total scalar dimension according * to those flat indices dim(V). - * 4.Each entry of RE is associated with a flat index k via FlatVectorBackend. + * 4.Each entry of RE is associated with a flat index k via flatVectorView. * This is normally a lexicographic index. The total scalar dimension according * to those flat indices dim(RE). * 5.Via those flat indices we now interpret C,V, and RE as vectors and compute the diadic @@ -117,9 +117,7 @@ public: template<typename Node, typename TreePath> void leaf(Node& node, TreePath treePath) { - using LocalBasisRange = typename Node::FiniteElement::Traits::LocalBasisType::Traits::RangeType; using MultiIndex = typename LocalIndexSet::MultiIndex; - using CoefficientBlock = typename std::decay<decltype(std::declval<Vector>()[std::declval<MultiIndex>()])>::type; using RangeBlock = typename std::decay<decltype(nodeToRangeEntry_(node, y_))>::type; auto&& fe = node.finiteElement(); @@ -129,7 +127,7 @@ public: localBasis.evaluateFunction(x_, shapeFunctionValues); // Get range entry associated to this node - auto&& re = nodeToRangeEntry_(node, y_); + auto re = flatVectorView(nodeToRangeEntry_(node, y_)); for (size_type i = 0; i < localBasis.size(); ++i) @@ -137,26 +135,25 @@ public: auto&& multiIndex = localIndexSet_.index(node.localIndex(i)); // Get coefficient associated to i-th shape function - auto&& c = coefficients_[multiIndex]; + auto c = flatVectorView(coefficients_[multiIndex]); // Get value of i-th shape function - auto&& v = shapeFunctionValues[i]; + auto v = flatVectorView(shapeFunctionValues[i]); + + // Notice that the range entry re, the coefficient c, and the shape functions // value v may all be scalar, vector, matrix, or general container valued. // The matching of their entries is done via the multistage procedure described // in the class documentation of DiscreteGlobalBasisFunction. - auto dimC = FlatVectorBackend<CoefficientBlock>::size(c); - auto dimV = FlatVectorBackend<LocalBasisRange>::size(v); - assert(dimC*dimV == FlatVectorBackend<RangeBlock>::size(re)); + auto&& dimC = c.size(); + auto dimV = v.size(); + assert(dimC*dimV == re.size()); for(size_type j=0; j<dimC; ++j) { - auto&& c_j = FlatVectorBackend<CoefficientBlock>::getEntry(c, j); + auto&& c_j = c[j]; for(size_type k=0; k<dimV; ++k) - { - auto&& v_k = FlatVectorBackend<LocalBasisRange>::getEntry(v, k); - FlatVectorBackend<RangeBlock>::getEntry(re, j*dimV + k) += c_j*v_k; - } + re[j*dimV + k] += c_j*v[k]; } } }