1 #ifndef EIGENVECTOR_HH
2 #define EIGENVECTOR_HH
3 
4 #ifdef HAVE_EIGEN
5 
6 #include <algorithm>
7 #include <iostream>
8 
9 #include <dune/common/densevector.hh>
10 #include <dune/common/ftraits.hh>
11 
12 #include <Eigen/Dense>
13 
14 namespace Dune
15 {
16 
17   namespace Fem
18   {
19     // forward declaration
20     template< class K > class EigenVector;
21   }
22 
23   // specialization of DenseMatVecTraits for EigenVector
24   template< class K >
25   struct DenseMatVecTraits< Fem::EigenVector< K > >
26   {
27     typedef Fem::EigenVector< K > derived_type;
28     typedef Eigen::Matrix< K, Eigen::Dynamic, 1 > container_type;
29     typedef K value_type;
30     typedef unsigned int size_type;
31   };
32 
33   template< class K >
34   struct FieldTraits< Fem::EigenVector< K > >
35   {
36     typedef typename FieldTraits< K >::field_type field_type;
37     typedef typename FieldTraits< K >::real_type real_type;
38   };
39 
40 
41   namespace Fem
42   {
43 
44     /** \brief An implementation of DenseVector which uses Eigen::Matrix<K, Eigen::Dynamic, 1> to provide the fields
45      *
46      * \tparam K is the field type (use float, double, complex, etc)
47      */
48     template< class K >
49     class EigenVector : public DenseVector< EigenVector< K > >
50     {
51       typedef EigenVector< K > ThisType;
52       typedef DenseVector< ThisType > BaseType;
53 
54     public:
55       typedef typename BaseType::size_type size_type;
56       typedef size_type SizeType;
57       typedef typename BaseType::value_type value_type;
58       typedef value_type FieldType;
59       typedef typename DenseMatVecTraits< ThisType >::container_type container_type;
60       typedef container_type DofStorageType;
61 
62       //! Constructor setting up a vector of a specified size
EigenVector(size_type size=0)63       explicit EigenVector( size_type size = 0 )
64       : data_( size )
65       {}
66 
67       //! Constructor setting up a vector initialized with a constant value
EigenVector(size_type size,const value_type & s)68       EigenVector( size_type size, const value_type& s )
69       : data_( size )
70       {
71         std::fill( data_.begin(), data_.end(), s );
72       }
73 
74       //! Copy constructor setting up a vector with the data of another one
75       template< class T >
EigenVector(const DenseVector<T> & v)76       EigenVector( const DenseVector< T >& v )
77       : data_( v.size() )
78       {
79         std::copy( v.begin(), v.end(), data_.begin() );
80       }
81 
82       //! Copy assignment operator
operator =(const EigenVector & other)83       EigenVector &operator=( const EigenVector& other )
84       {
85         data_.resize( other.size() );
86         std::copy( other.begin(), other.end(), data_.begin() );
87         return *this;
88       }
89 
operator [](size_type index) const90       const value_type& operator[]( size_type index ) const
91       {
92         return data_( index );
93       }
94 
operator [](size_type index)95       value_type& operator[]( size_type index )
96       {
97         return data_( index );
98       }
99 
100       //! Obtain coefficients
coefficients() const101       const DofStorageType& coefficients() const
102       {
103         return data_;
104       }
105 
106       //! Obtain coefficients
coefficients()107       DofStorageType& coefficients()
108       {
109         return data_;
110       }
111 
112       //! Obtain pointer to data
data() const113       const value_type* data() const
114       {
115         return data_.data();
116       }
117 
118       //! Obtain pointer to data
data()119       value_type* data()
120       {
121         return data_.data();
122       }
123 
124       //! Allocate memory
reserve(size_type newSize)125       void reserve( size_type newSize )
126       {
127         data_.resize( newSize );
128       }
129 
130       //! Resize vector
resize(size_type newSize)131       void resize( size_type newSize )
132       {
133         data_.resize( newSize );
134       }
135 
136       //! Resize vector and initialize
resize(size_type newSize,const value_type & s)137       void resize( size_type newSize, const value_type& s )
138       {
139         data_.resize( newSize );
140         std::fill( data_.begin(), data_.end(), s );
141       }
142 
size() const143       size_type size() const
144       {
145         return data_.size();
146       }
147 
148     private:
149       DofStorageType data_;
150     };
151 
152     /** \brief Read a EigenVector from an input stream
153      *  \relates EigenVector
154      *
155      *  \note This operator is STL compilant, i.e., the content of v is only
156      *        changed if the read operation is successful.
157      *
158      *  \param[in]  in  std::istream to read from
159      *  \param[out] v   EigenVector to be read
160      *
161      *  \returns the input stream (in)
162      */
163     template< class K >
operator >>(std::istream & in,EigenVector<K> & v)164     inline std::istream& operator>>( std::istream& in, EigenVector< K >& v )
165     {
166       EigenVector< K > w(v);
167       for( typename EigenVector< K >::size_type i = 0; i < w.size(); ++i )
168         in >> w[ i ];
169       if(in)
170         v = w;
171       return in;
172     }
173 
174   } // namespace Fem
175 
176 } // namespace Dune
177 
178 #endif
179 
180 #endif // #ifndef EIGENVECTOR_HH
181