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