dune-pdelab  2.4.1
eigen/vector.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
2 #define DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
3 
4 #if HAVE_EIGEN
5 
6 #include <memory>
7 
8 #include <dune/istl/bvector.hh>
9 
13 #include "descriptors.hh"
14 #include <Eigen/Dense>
15 
16 namespace Dune {
17  namespace PDELab {
18 
19  namespace EIGEN {
20 
21  template<typename GFS, typename ET>
22  class VectorContainer
23  : public Backend::impl::Wrapper<Eigen::Matrix<ET, Eigen::Dynamic, 1>>
24  {
25  public:
26  typedef Eigen::Matrix<ET, Eigen::Dynamic, 1> Container;
27 
28  private:
29 
30  friend Backend::impl::Wrapper<Container>;
31 
32  public:
33  typedef ET ElementType;
34  typedef ET E;
35 
36  // for ISTL solver compatibility
37  typedef ElementType field_type;
38 
39  typedef GFS GridFunctionSpace;
40  typedef std::size_t size_type;
41 
42  typedef typename GFS::Ordering::Traits::ContainerIndex ContainerIndex;
43 
44  typedef ElementType* iterator;
45  typedef const ElementType* const_iterator;
46  // #warning iterators does not work well with Eigen
47  // typedef typename Container::iterator iterator;
48  // typedef typename Container::const_iterator const_iterator;
49 
50  template<typename LFSCache>
51  using LocalView = UncachedVectorView<VectorContainer,LFSCache>;
52 
53  template<typename LFSCache>
54  using ConstLocalView = ConstUncachedVectorView<const VectorContainer,LFSCache>;
55 
56  VectorContainer(const VectorContainer& rhs)
57  : _gfs(rhs._gfs)
58  , _container(std::make_shared<Container>(*(rhs._container)))
59  {}
60 
61  VectorContainer (const GFS& gfs, Backend::attached_container = Backend::attached_container())
62  : _gfs(gfs)
63  , _container(std::make_shared<Container>(gfs.ordering().blockCount()))
64  {}
65 
67  VectorContainer(const GFS& gfs, Backend::unattached_container)
68  : _gfs(gfs)
69  {}
70 
76  VectorContainer (const GFS& gfs, Container& container)
77  : _gfs(gfs)
78  , _container(stackobject_to_shared_ptr(container))
79  {
80  _container->resize(gfs.ordering().blockCount());
81  }
82 
83  VectorContainer (const GFS& gfs, const E& e)
84  : _gfs(gfs)
85  , _container(std::make_shared<Container>(Container::Constant(gfs.ordering().blockCount(),e)))
86  {}
87 
88  void detach()
89  {
90  _container.reset();
91  }
92 
93  void attach(std::shared_ptr<Container> container)
94  {
95  _container = container;
96  }
97 
98  bool attached() const
99  {
100  return bool(_container);
101  }
102 
103  const std::shared_ptr<Container>& storage() const
104  {
105  return _container;
106  }
107 
108  size_type N() const
109  {
110  return _container->size();
111  }
112 
113  VectorContainer& operator=(const VectorContainer& r)
114  {
115  if (this == &r)
116  return *this;
117  if (attached())
118  {
119  (*_container) = (*r._container);
120  }
121  else
122  {
123  _container = std::make_shared<Container>(*(r._container));
124  }
125  return *this;
126  }
127 
128  VectorContainer& operator=(const E& e)
129  {
130  (*_container) = Container::Constant(N(),e);
131  return *this;
132  }
133 
134  VectorContainer& operator*=(const E& e)
135  {
136  (*_container) *= e;
137  return *this;
138  }
139 
140 
141  VectorContainer& operator+=(const E& e)
142  {
143  (*_container) += Container::Constant(N(),e);
144  return *this;
145  }
146 
147  VectorContainer& operator+=(const VectorContainer& y)
148  {
149  (*_container) += (*y._container);
150  return *this;
151  }
152 
153  VectorContainer& operator-= (const VectorContainer& y)
154  {
155  (*_container) -= (*y._container);
156  return *this;
157  }
158 
159  E& operator[](const ContainerIndex& ci)
160  {
161  return (*_container)(ci[0]);
162  }
163 
164  const E& operator[](const ContainerIndex& ci) const
165  {
166  return (*_container)(ci[0]);
167  }
168 
169  typename Dune::template FieldTraits<E>::real_type two_norm() const
170  {
171  return _container->norm();
172  }
173 
174  typename Dune::template FieldTraits<E>::real_type one_norm() const
175  {
176  return _container->template lpNorm<1>();
177  }
178 
179  typename Dune::template FieldTraits<E>::real_type infinity_norm() const
180  {
181  return _container->template lpNorm<Eigen::Infinity>();
182  }
183 
185  E operator*(const VectorContainer& y) const
186  {
187  return _container->transpose() * (*y._container);
188  }
189 
191  E dot(const VectorContainer& y) const
192  {
193  return _container->dot(*y._container);
194  }
195 
197  VectorContainer& axpy(const E& a, const VectorContainer& y)
198  {
199  (*_container) += a * (*y._container);
200  return *this;
201  }
202 
203  // for debugging and AMG access
204  Container& base ()
205  {
206  return *_container;
207  }
208 
209  const Container& base () const
210  {
211  return *_container;
212  }
213 
214  private:
215 
216  Container& native ()
217  {
218  return *_container;
219  }
220 
221  const Container& native () const
222  {
223  return *_container;
224  }
225 
226  public:
227 
228  iterator begin()
229  {
230  return _container->data();
231  }
232 
233  const_iterator begin() const
234  {
235  return _container->data();
236  }
237 
238  iterator end()
239  {
240  return _container->data() + N();
241  }
242 
243  const_iterator end() const
244  {
245  return _container->data() + N();
246  }
247 
248  size_t flatsize() const
249  {
250  return _container->size();
251  }
252 
253  const GFS& gridFunctionSpace() const
254  {
255  return _gfs;
256  }
257 
258  private:
259  const GFS& _gfs;
260  std::shared_ptr<Container> _container;
261 
262  };
263 
264  } // end namespace EIGEN
265 
266 
267 #ifndef DOXYGEN
268 
269  template<typename GFS, typename E>
270  struct EigenVectorSelectorHelper
271  {
272  using Type = EIGEN::VectorContainer<GFS, E>;
273  };
274 
275  namespace Backend {
276  namespace impl {
277 
278  template<typename GFS, typename E>
279  struct BackendVectorSelectorHelper<EigenVectorBackend, GFS, E>
280  : public EigenVectorSelectorHelper<GFS,E>
281  {};
282 
283  } // namespace impl
284  } // namespace Backend
285 
286 #endif // DOXYGEN
287 
288  } // namespace PDELab
289 } // namespace Dune
290 
291 #endif
292 
293 #endif // DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
294 
295 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
296 // vi: set et ts=4 sw=2 sts=2:
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:199
STL namespace.
Definition: adaptivity.hh:27
Backend::unattached_container unattached_container
Definition: backend/common/tags.hh:38
Backend::attached_container attached_container
Definition: backend/common/tags.hh:39
const Entity & e
Definition: localfunctionspace.hh:111
Various tags for influencing backend behavior.