3#ifndef DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
4#define DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
14#include <dune/pdelab/backend/interface.hh>
15#include <dune/pdelab/backend/common/uncachedvectorview.hh>
16#include "descriptors.hh"
29 template<
typename GFS,
typename ET>
31 :
public Backend::impl::Wrapper<::Eigen::Matrix<ET, ::Eigen::Dynamic, 1>>
34 typedef ::Eigen::Matrix<ET, ::Eigen::Dynamic, 1> Container;
38 friend Backend::impl::Wrapper<Container>;
41 typedef ET ElementType;
45 typedef ElementType field_type;
47 typedef GFS GridFunctionSpace;
48 typedef std::size_t size_type;
50 typedef typename GFS::Ordering::Traits::ContainerIndex ContainerIndex;
52 typedef ElementType* iterator;
53 typedef const ElementType* const_iterator;
58 template<
typename LFSCache>
59 using LocalView = UncachedVectorView<VectorContainer,LFSCache>;
61 template<
typename LFSCache>
62 using ConstLocalView = ConstUncachedVectorView<const VectorContainer,LFSCache>;
64 VectorContainer(
const VectorContainer& rhs)
66 , _container(
std::make_shared<Container>(*(rhs._container)))
69 VectorContainer (
const GFS& gfs, Backend::attached_container = Backend::attached_container())
71 , _container(
std::make_shared<Container>(gfs.ordering().blockCount()))
75 VectorContainer(
const GFS& gfs, Backend::unattached_container)
84 VectorContainer (
const GFS& gfs, Container& container)
88 _container->resize(gfs.ordering().blockCount());
91 VectorContainer (
const GFS& gfs,
const E& e)
93 , _container(
std::make_shared<Container>(Container::Constant(gfs.ordering().blockCount(),e)))
101 void attach(std::shared_ptr<Container> container)
103 _container = container;
106 bool attached()
const
108 return bool(_container);
111 const std::shared_ptr<Container>& storage()
const
118 return _container->size();
121 VectorContainer& operator=(
const VectorContainer& r)
127 (*_container) = (*r._container);
131 _container = std::make_shared<Container>(*(r._container));
136 VectorContainer& operator=(
const E& e)
138 (*_container) = Container::Constant(N(),e);
142 VectorContainer& operator*=(
const E& e)
149 VectorContainer& operator+=(
const E& e)
151 (*_container) += Container::Constant(N(),e);
155 VectorContainer& operator+=(
const VectorContainer& y)
157 (*_container) += (*y._container);
161 VectorContainer& operator-= (
const VectorContainer& y)
163 (*_container) -= (*y._container);
167 E& operator[](
const ContainerIndex& ci)
169 return (*_container)(ci[0]);
172 const E& operator[](
const ContainerIndex& ci)
const
174 return (*_container)(ci[0]);
177 typename Dune::template FieldTraits<E>::real_type two_norm()
const
179 return _container->norm();
182 typename Dune::template FieldTraits<E>::real_type one_norm()
const
184 return _container->template lpNorm<1>();
187 typename Dune::template FieldTraits<E>::real_type infinity_norm()
const
189 return _container->template lpNorm<::Eigen::Infinity>();
193 E operator*(
const VectorContainer& y)
const
195 return _container->transpose() * (*y._container);
199 E
dot(
const VectorContainer& y)
const
201 return _container->dot(*y._container);
205 VectorContainer& axpy(
const E& a,
const VectorContainer& y)
207 (*_container) += a * (*y._container);
217 const Container& base ()
const
229 const Container& native ()
const
238 return _container->data();
241 const_iterator begin()
const
243 return _container->data();
248 return _container->data() + N();
251 const_iterator end()
const
253 return _container->data() + N();
256 size_t flatsize()
const
258 return _container->size();
261 const GFS& gridFunctionSpace()
const
268 std::shared_ptr<Container> _container;
277 template<
typename GFS,
typename E>
278 struct EigenVectorSelectorHelper
280 using Type = PDELab::Eigen::VectorContainer<GFS, E>;
286 template<
typename GFS,
typename E>
287 struct BackendVectorSelectorHelper<Eigen::VectorBackend, GFS, E>
288 :
public EigenVectorSelectorHelper<GFS,E>
This file implements a vector space as a tensor product of a given vector space. The number of compon...
auto dot(const A &a, const B &b) -> typename std::enable_if<!IsVector< A >::value &&!std::is_same< typename FieldTraits< A >::field_type, typename FieldTraits< A >::real_type > ::value, decltype(conj(a) *b)>::type
computes the dot product for fundamental data types according to Petsc's VectDot function: dot(a,...
Definition: dotproduct.hh:40
Dune namespace.
Definition: alignedallocator.hh:11
std::shared_ptr< T > stackobject_to_shared_ptr(T &t)
Create a shared_ptr for a stack-allocated object.
Definition: shared_ptr.hh:70
This file implements several utilities related to std::shared_ptr.