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