dune-fem  2.8-git
space/common/dataprojection/default.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_HPDG_SPACE_COMMON_DATAPROJECTION_DEFAULT_HH
2 #define DUNE_FEM_HPDG_SPACE_COMMON_DATAPROJECTION_DEFAULT_HH
3 
4 #include <cassert>
5 #include <cstddef>
6 
7 #include <functional>
8 #include <vector>
9 
10 #include <dune/common/dynvector.hh>
11 
16 
17 #include "dataprojection.hh"
18 
19 namespace Dune
20 {
21 
22  namespace Fem
23  {
24 
25  namespace hpDG
26  {
27 
28  // DefaultDataProjection
29  // ---------------------
30 
37  template< class DiscreteFunction >
39  : public DataProjection< typename DiscreteFunction::DiscreteFunctionSpaceType, DefaultDataProjection< DiscreteFunction > >
40  {
43 
44  public:
50  using EntityType = typename BaseType::EntityType;
51 
52  private:
53  using RangeFieldType = typename DiscreteFunction::RangeFieldType;
54  using LocalDofVectorType = Dune::DynamicVector< RangeFieldType >;
56 
57  static const std::size_t localBlockSize = DiscreteFunctionSpaceType::localBlockSize;
59 
60  public:
65  explicit DefaultDataProjection ( DiscreteFunction &discreteFunction )
66  : discreteFunction_( discreteFunction )
67  {
68  discreteFunction_.get().enableDofCompression();
69  }
70 
73 #ifndef DOXYGEN
74 
75  DefaultDataProjection ( const ThisType & ) = delete;
76 
77  DefaultDataProjection ( ThisType && ) = default;
78 
79  ThisType &operator= ( const ThisType & ) = delete;
80 
81  ThisType &operator= ( ThisType && ) = default;
82 
83 #endif // #ifndef DOXYGEN
84 
86  void operator() ( const EntityType &entity,
87  const BasisFunctionSetType &prior,
88  const BasisFunctionSetType &present,
89  const std::vector< std::size_t > &origin,
90  const std::vector< std::size_t > &destination )
91  {
92  LocalFunctionType localFunc( prior );
93  read( origin, localFunc.localDofVector() );
94 
95  assert( present.size() == space().basisFunctionSet( entity ).size() );
96  LocalDofVectorType localDofVector( present.size() );
97 
98  LocalInterpolationType& ip = this->interpolation();
99  auto guard = bindGuard( ip, entity );
100  ip( localFunc, localDofVector );
101 
102  write( destination, localDofVector );
103  }
104 
106  template <class TemporaryStorage>
107  void operator () ( TemporaryStorage& tmp )
108  {
109  auto& df = discreteFunction();
110 
111  // copy dofs to temporary storage, old order still exists
112  auto dfit = df.dbegin();
113  const auto endtmp = tmp.dend();
114  for( auto it = tmp.dbegin(); it != endtmp; ++it, ++dfit )
115  {
116  assert( dfit != df.dend() );
117  *it = *dfit;
118  }
119 
122  LocalInterpolationType interpolation( df.space() );
123 
124  // interpolate to new space, this can be a
125  // Lagrange interpolation or a L2 projection, both locally
126  const auto end = df.space().end();
127  for( auto it = df.space().begin(); it != end; ++it )
128  {
129  const auto& entity = *it;
130 
131  auto tguard = bindGuard( tmpLF, entity );
132  auto lguard = bindGuard( lf, entity );
133 
134  // if size is the same we can just copy the dof values
135  if( tmpLF.size() == lf.size() )
136  {
137  lf.assign( tmpLF );
138  }
139  else
140  {
141  auto iguard = bindGuard( interpolation, entity );
142  // otherwise a local interpolation is needed
143  interpolation( tmpLF, lf );
144  }
145  }
146  }
147 
149  template< class Communicator >
150  void addToList ( Communicator &comm )
151  {
152  comm.addToList( discreteFunction() );
153  }
154 
155  private:
156  template< class LocalDofVector >
157  void read ( const std::vector< std::size_t > &blocks, LocalDofVector &localDofVector ) const
158  {
159  assert( localDofVector.size() == localBlockSize*blocks.size() );
160  std::size_t index = 0;
161  for( auto i : blocks )
162  {
163  const auto block = discreteFunction().block( i );
164  for( std::size_t j = 0; j < localBlockSize; ++j )
165  localDofVector[ index++ ] = (*block)[ j ];
166  }
167  assert( index == localDofVector.size() );
168  }
169 
170  template< class LocalDofVector >
171  void write ( const std::vector< std::size_t > &blocks, const LocalDofVector &localDofVector )
172  {
173  assert( localDofVector.size() == localBlockSize*blocks.size() );
174  std::size_t index = 0;
175  for( auto i : blocks )
176  {
177  auto block = discreteFunction().block( i );
178  for( std::size_t j = 0; j < localBlockSize; ++j )
179  (*block)[ j ] = localDofVector[ index++ ];
180  }
181  assert( index == localDofVector.size() );
182  }
183 
184  DiscreteFunction &discreteFunction () { return discreteFunction_.get(); }
185 
186  const DiscreteFunction &discreteFunction () const { return discreteFunction_.get(); }
187 
188  const DiscreteFunctionSpaceType &space () const { return discreteFunction().space(); }
189  LocalInterpolationType& interpolation ()
190  {
191  if( ! interpolation_ )
192  interpolation_.reset( new LocalInterpolationType( space() ) );
193  return *interpolation_;
194  }
195 
196 
197  std::reference_wrapper< DiscreteFunction > discreteFunction_;
198  std::shared_ptr< LocalInterpolationType > interpolation_;
199  };
200 
201  } // namespace hpDG
202 
203  // forward types to Fem namespace for convenience
204  using hpDG::DefaultDataProjection;
205 
206  } // namespace Fem
207 
208 } // namespace Dune
209 
210 #endif // #ifndef DUNE_FEM_HPDG_SPACE_COMMON_DATAPROJECTION_DEFAULT_HH
Definition: bindguard.hh:11
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:517
static auto bindGuard(Object &object, Args &&... args) -> std::enable_if_t< isBindable< Object, Args... >::value, BindGuard< Object > >
Definition: bindguard.hh:67
interface for local functions
Definition: localfunction.hh:77
const LocalDofVectorType & localDofVector() const
return const reference to local Dof Vector
Definition: localfunction.hh:415
void assign(const LocalFunction< BasisFunctionSet, T > &other)
assign all DoFs of this local function
Definition: localfunction.hh:189
SizeType size() const
obtain the number of local DoFs
Definition: localfunction.hh:360
Definition: mutable.hh:31
Abstract definition of the local restriction and prolongation of discrete functions.
Definition: dataprojection/dataprojection.hh:29
typename BasisFunctionSetType::EntityType EntityType
entity type
Definition: dataprojection/dataprojection.hh:38
DiscreteFunctionSpace DiscreteFunctionSpaceType
discrete function space type
Definition: dataprojection/dataprojection.hh:34
typename DiscreteFunctionSpaceType::BasisFunctionSetType BasisFunctionSetType
basis function set type
Definition: dataprojection/dataprojection.hh:36
Local -projection for the restriction and prolongation of discrete functions.
Definition: space/common/dataprojection/default.hh:40
void addToList(Communicator &comm)
()
Definition: space/common/dataprojection/default.hh:150
typename BaseType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
discrete function space type
Definition: space/common/dataprojection/default.hh:46
typename BaseType::EntityType EntityType
entity type
Definition: space/common/dataprojection/default.hh:50
typename BaseType::BasisFunctionSetType BasisFunctionSetType
basis function set type
Definition: space/common/dataprojection/default.hh:48
void operator()(const EntityType &entity, const BasisFunctionSetType &prior, const BasisFunctionSetType &present, const std::vector< std::size_t > &origin, const std::vector< std::size_t > &destination)
Definition: space/common/dataprojection/default.hh:86
DefaultDataProjection(DiscreteFunction &discreteFunction)
Definition: space/common/dataprojection/default.hh:65
Definition: common/localinterpolation.hh:22