3 #ifndef DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
4 #define DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
10 #include <dune/istl/bvector.hh>
15 #include "descriptors.hh"
16 #include <Eigen/Dense>
28 template<
typename GFS,
typename ET>
30 :
public Backend::impl::Wrapper<::Eigen::Matrix<ET, ::Eigen::Dynamic, 1>>
33 typedef ::Eigen::Matrix<ET, ::Eigen::Dynamic, 1> Container;
37 friend Backend::impl::Wrapper<Container>;
40 typedef ET ElementType;
44 typedef ElementType field_type;
46 typedef GFS GridFunctionSpace;
47 typedef std::size_t size_type;
49 typedef typename GFS::Ordering::Traits::ContainerIndex ContainerIndex;
51 typedef ElementType* iterator;
52 typedef const ElementType* const_iterator;
57 template<
typename LFSCache>
58 using LocalView = UncachedVectorView<VectorContainer,LFSCache>;
60 template<
typename LFSCache>
61 using ConstLocalView = ConstUncachedVectorView<const VectorContainer,LFSCache>;
63 VectorContainer(
const VectorContainer& rhs)
65 , _container(std::make_shared<Container>(*(rhs._container)))
68 VectorContainer (
const GFS& gfs, Backend::attached_container = Backend::attached_container())
70 , _container(std::make_shared<Container>(gfs.ordering().blockCount()))
74 VectorContainer(
const GFS& gfs, Backend::unattached_container)
83 VectorContainer (
const GFS& gfs, Container& container)
85 , _container(stackobject_to_shared_ptr(container))
87 _container->resize(gfs.ordering().blockCount());
90 VectorContainer (
const GFS& gfs,
const E&
e)
92 , _container(std::make_shared<Container>(Container::Constant(gfs.ordering().blockCount(),
e)))
100 void attach(std::shared_ptr<Container> container)
102 _container = container;
105 bool attached()
const
107 return bool(_container);
110 const std::shared_ptr<Container>& storage()
const
117 return _container->size();
120 VectorContainer& operator=(
const VectorContainer& r)
126 (*_container) = (*r._container);
130 _container = std::make_shared<Container>(*(r._container));
135 VectorContainer& operator=(
const E&
e)
137 (*_container) = Container::Constant(N(),
e);
141 VectorContainer& operator*=(
const E&
e)
148 VectorContainer& operator+=(
const E&
e)
150 (*_container) += Container::Constant(N(),
e);
154 VectorContainer& operator+=(
const VectorContainer& y)
156 (*_container) += (*y._container);
160 VectorContainer& operator-= (
const VectorContainer& y)
162 (*_container) -= (*y._container);
166 E& operator[](
const ContainerIndex& ci)
168 return (*_container)(ci[0]);
171 const E& operator[](
const ContainerIndex& ci)
const
173 return (*_container)(ci[0]);
176 typename Dune::template FieldTraits<E>::real_type two_norm()
const
178 return _container->norm();
181 typename Dune::template FieldTraits<E>::real_type one_norm()
const
183 return _container->template lpNorm<1>();
186 typename Dune::template FieldTraits<E>::real_type infinity_norm()
const
188 return _container->template lpNorm<::Eigen::Infinity>();
192 E operator*(
const VectorContainer& y)
const
194 return _container->transpose() * (*y._container);
198 E dot(
const VectorContainer& y)
const
200 return _container->dot(*y._container);
204 VectorContainer& axpy(
const E& a,
const VectorContainer& y)
206 (*_container) += a * (*y._container);
216 const Container& base ()
const
228 const Container&
native ()
const
237 return _container->data();
240 const_iterator begin()
const
242 return _container->data();
247 return _container->data() + N();
250 const_iterator end()
const
252 return _container->data() + N();
255 size_t flatsize()
const
257 return _container->size();
260 const GFS& gridFunctionSpace()
const
267 std::shared_ptr<Container> _container;
276 template<
typename GFS,
typename E>
277 struct EigenVectorSelectorHelper
279 using Type = PDELab::Eigen::VectorContainer<GFS, E>;
285 template<
typename GFS,
typename E>
286 struct BackendVectorSelectorHelper<Eigen::VectorBackend, GFS, E>
287 :
public EigenVectorSelectorHelper<GFS,E>
300 #endif // DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH