constantspace.hh

Go to the documentation of this file.

/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
/*                                                                           */
/*  This file is part of the library KASKADE 7                               */
/*  https://www.zib.de/research/projects/kaskade7-finite-element-toolbox     */
/*                                                                           */
/*  Copyright (C) 2002-2023 Zuse Institute Berlin                            */
/*                                                                           */
/*  KASKADE 7 is distributed under the terms of the ZIB Academic License.    */
/*    see $KASKADE/academic.txt                                              */
/*                                                                           */
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
 
#ifndef CONSTANTSPACE_HH
#define CONSTANTSPACE_HH
 
#include <numeric>
#include <limits>
 
#include "fem/lagrangeshapefunctions.hh"
#include "fem/functionspace.hh"
#include "fem/views.hh"
 
namespace Kaskade
{
  template <class ScalarType, class GV>
  class ConstantMapper
  {
    typedef ConstantMapper<ScalarType,GV> Self;
 
  public:
    typedef ScalarType Scalar;
 
    typedef GV GridView;
    typedef typename GridView::Grid Grid;
    typedef typename GridView::IndexSet IndexSet;
    using ctype = typename GridView::ctype;
 
    static int const dim = Grid::dimension;
    static int const maxLocalSize = 1;
 
    using Cell = Kaskade::Cell<GridView>;
    typedef LagrangeSimplexShapeFunctionSet<ctype,dim,ScalarType> ShapeFunctionSetImplementation;  // constant functions are the same for any polytope
    typedef typename LagrangeShapeFunctionSetContainer<ctype,dim,ScalarType>::value_type ShapeFunctionSet;
    typedef RangeView<std::vector<size_t>::const_iterator> GlobalIndexRange;
 
    typedef std::pair<size_t,int>                  IndexPair;
    typedef std::vector<IndexPair>::const_iterator SortedIndexIterator;
    typedef RangeView<SortedIndexIterator>         SortedIndexRange;
  
    static int const continuity = std::numeric_limits<int>::max()-1; // mpl::int_ chokes on true max int value...
    
    static bool const globalSupport = true;
 
    template <int m>
    struct Element
    {
      typedef FunctionSpaceElement<FEFunctionSpace<Self>,m> type;
    };
 
    ConstantMapper(GridView const& gridView)
    : indexSet(gridView.indexSet())
    , globIdx{0}
    , sortedIdx(1)
    {
      sortedIdx[0] = std::make_pair(0,0);
    }
 
    virtual ~ConstantMapper() {}
 
    ShapeFunctionSet const& shapefunctions(Cell const& cell) const
    {
      return lagrangeShapeFunctionSet<ctype,dim,Scalar>(Dune::GeometryType(Dune::GeometryType::simplex,dim),0);  // constant functions are the same for any polytope
    }
 
    ShapeFunctionSet const& shapefunctions(size_t) const
    {
      // constant functions are the same for any polytope
      Dune::GeometryType gt(0,dim,false);     // elementId:0, dimension:dim, bool isNone:false
      return lagrangeShapeFunctionSet<ctype,dim,Scalar>(gt,0);
    }
 
    int maxOrder() const
    {
      return 0;
    }
 
    GlobalIndexRange initGlobalIndexRange() const 
    {
      return GlobalIndexRange(globIdx.begin(),globIdx.begin());
    }
 
    GlobalIndexRange globalIndices(typename Grid::template Codim<0>::Entity const& /* cell */) const
    {
      return GlobalIndexRange(globIdx.begin(),globIdx.end());
    }
 
    GlobalIndexRange globalIndices(size_t) const 
    {
      return GlobalIndexRange(globIdx.begin(),globIdx.end());
    }
    
    static SortedIndexRange initSortedIndexRange() 
    {
      static std::vector<IndexPair> dummy; // empty
      return SortedIndexRange(dummy.begin(),dummy.end());
    }
    
    SortedIndexRange sortedIndices(Cell const& ) const 
    {
      return SortedIndexRange(sortedIdx.begin(),sortedIdx.end());
    }
    
    SortedIndexRange sortedIndices(size_t) const 
    {
      return SortedIndexRange(sortedIdx.begin(),sortedIdx.end());
    }
 
    size_t size() const { return 1; }
 
    class Converter
    {
    public:
      Converter() {}
 
      Converter(Cell const& /* cell */) {}
 
      void moveTo(Cell const& /* cell */) { }
 
      void setLocalPosition(Dune::FieldVector<ctype,Grid::dimension> const& /* x */) {}
 
      Dune::FieldMatrix<Scalar,1,1> global(Dune::FieldMatrix<Scalar,1,1> const& sf) const { return sf; }
 
 
      template <class Scalar>
      VariationalArg<Scalar,dim,1> global(VariationalArg<Scalar,dim,1>  const& sf, int deriv=1) const
      {
        VariationalArg<Scalar,dim,1> phi;
        phi.value = 1;
        phi.derivative = 0;
        phi.hessian = 0;
        return phi;
      }
        
      Dune::FieldMatrix<Scalar,1,1> local(Dune::FieldMatrix<Scalar,1,1> const& glob) const { return glob; }
    };
 
    class Combiner {
    public:
      Combiner(Cell const& cell) {}
 
      template <class Matrix>
      void rightTransform(Matrix& /* A */) const {}
 
      template <int n, int m>
      void rightTransform(std::vector<VariationalArg<Scalar,n,m> >& /* v */) const {}
 
      template <class Matrix>
      void leftPseudoInverse(Matrix& /* A */) const {}
 
      operator Dune::BCRSMatrix<Dune::FieldMatrix<Scalar,1,1> >() const
      {
        Dune::BCRSMatrix<Dune::FieldMatrix<Scalar,1,1> > K(1,1,Dune::BCRSMatrix<Dune::FieldMatrix<Scalar,1,1> >::random);
        K.incrementrowsize(0);
        K.endrowsizes();
        K.addindex(0,0);
        K.endindices();
 
        *K[0].begin() = 1;
        return K;
      }
 
    };
 
    Combiner combiner(Cell const& cell, size_t index) const { return Combiner(cell); }
 
 
 
    void update() {}
 
  private:
    IndexSet const& indexSet;
    std::vector<size_t> globIdx;
    std::vector<IndexPair> sortedIdx;
  };
} /* end of namespace Kaskade */
 
#endif