diff --git a/dune/gfe/parallel/p2mapper.hh b/dune/gfe/parallel/p2mapper.hh new file mode 100644 index 0000000000000000000000000000000000000000..d6eed91ebe341a57b976ac1e2707cb2f09932e4e --- /dev/null +++ b/dune/gfe/parallel/p2mapper.hh @@ -0,0 +1,57 @@ +// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- +// vi: set et ts=4 sw=2 sts=2: +#ifndef DUNE_GFE_PARALLEL_P2MAPPER_HH +#define DUNE_GFE_PARALLEL_P2MAPPER_HH + +#include <dune/geometry/type.hh> +#include <dune/common/typetraits.hh> + +#include <dune/functions/functionspacebases/pqknodalbasis.hh> + +/** \brief Mimic a dune-grid mapper for a P2 space, using the dune-functions dof ordering of such a space + */ +template<class GridView> +class P2BasisMapper +{ + typedef typename GridView::Grid::template Codim<0>::Entity Element; +public: + + typedef typename Dune::Functions::PQkNodalBasis<GridView,2>::MultiIndex::value_type Index; + + P2BasisMapper(const GridView& gridView) + : p2Basis_(gridView) + {} + + std::size_t size() const + { + return p2Basis_.size(); + } + + template <class Entity> + bool contains(const Entity& entity, uint subEntity, uint codim, Index& result) const + { + auto localView = p2Basis_.localView(); + auto localIndexSet = p2Basis_.localIndexSet(); + localView.bind(entity); + localIndexSet.bind(localView); + + Index localIndex; + for (size_t i=0; i<localIndexSet.size(); i++) + { + if (localView.tree().finiteElement().localCoefficients().localKey(i).subEntity() == subEntity + and localView.tree().finiteElement().localCoefficients().localKey(i).codim() == codim) + { + localIndex = i; + result = localIndexSet.index(i)[0]; + return true; + } + } + + return false; + } + + Dune::Functions::PQkNodalBasis<GridView,2> p2Basis_; +}; + +#endif + diff --git a/dune/gfe/riemanniantrsolver.hh b/dune/gfe/riemanniantrsolver.hh index a19cda9346a39d060e764f394a5baadaf81f7b43..8a5c6f48ae4d4084dfa5c76c3fb670f496b6a34f 100644 --- a/dune/gfe/riemanniantrsolver.hh +++ b/dune/gfe/riemanniantrsolver.hh @@ -13,29 +13,27 @@ #include <dune/solvers/solvers/iterativesolver.hh> #include <dune/solvers/solvers/loopsolver.hh> -#include <dune/fufem/functionspacebases/p2nodalbasis.hh> - #include "geodesicfeassembler.hh" #include <dune/grid/utility/globalindexset.hh> #include <dune/gfe/parallel/globalp1mapper.hh> #include <dune/gfe/parallel/globalp2mapper.hh> +#include <dune/gfe/parallel/p2mapper.hh> /** \brief Assign GlobalMapper and LocalMapper types to a dune-fufem FunctionSpaceBasis */ -template <typename Basis> +template <typename GridView, typename Basis> struct MapperFactory {}; /** \brief Specialization for PQ1NodalBasis */ template <typename GridView> -struct MapperFactory<DuneFunctionsBasis<Dune::Functions::PQkNodalBasis<GridView,1> > > +struct MapperFactory<GridView, Dune::Functions::PQkNodalBasis<GridView,1> > { typedef Dune::GlobalP1Mapper<GridView> GlobalMapper; typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGVertexLayout> LocalMapper; }; -/** \brief Specialization for PQ2NodalBasis */ template <typename GridView> -struct MapperFactory<DuneFunctionsBasis<Dune::Functions::PQkNodalBasis<GridView,2> > > +struct MapperFactory<GridView, Dune::Functions::PQkNodalBasis<GridView,2> > { typedef Dune::GlobalP2Mapper<GridView> GlobalMapper; typedef P2BasisMapper<GridView> LocalMapper; @@ -43,7 +41,7 @@ struct MapperFactory<DuneFunctionsBasis<Dune::Functions::PQkNodalBasis<GridView, /** \brief Specialization for PQ3NodalBasis */ template <typename GridView> -struct MapperFactory<DuneFunctionsBasis<Dune::Functions::PQkNodalBasis<GridView,3> > > +struct MapperFactory<GridView, Dune::Functions::PQkNodalBasis<GridView,3> > { // Error: we don't currently have a global P3 mapper }; @@ -69,8 +67,8 @@ class RiemannianTrustRegionSolver typedef std::vector<TargetSpace> SolutionType; #if HAVE_MPI - typedef typename MapperFactory<Basis>::GlobalMapper GlobalMapper; - typedef typename MapperFactory<Basis>::LocalMapper LocalMapper; + typedef typename MapperFactory<typename Basis::GridView,Basis>::GlobalMapper GlobalMapper; + typedef typename MapperFactory<typename Basis::GridView,Basis>::LocalMapper LocalMapper; #endif