DenseStorage.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr>
5 // Copyright (C) 2006-2009 Benoit Jacob <jacob.benoit.1@gmail.com>
6 // Copyright (C) 2010 Hauke Heibel <hauke.heibel@gmail.com>
7 //
8 // Eigen is free software; you can redistribute it and/or
9 // modify it under the terms of the GNU Lesser General Public
10 // License as published by the Free Software Foundation; either
11 // version 3 of the License, or (at your option) any later version.
12 //
13 // Alternatively, you can redistribute it and/or
14 // modify it under the terms of the GNU General Public License as
15 // published by the Free Software Foundation; either version 2 of
16 // the License, or (at your option) any later version.
17 //
18 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
19 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
21 // GNU General Public License for more details.
22 //
23 // You should have received a copy of the GNU Lesser General Public
24 // License and a copy of the GNU General Public License along with
25 // Eigen. If not, see <http://www.gnu.org/licenses/>.
26 
27 #ifndef EIGEN_MATRIXSTORAGE_H
28 #define EIGEN_MATRIXSTORAGE_H
29 
30 #ifdef EIGEN_DENSE_STORAGE_CTOR_PLUGIN
31  #define EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN EIGEN_DENSE_STORAGE_CTOR_PLUGIN;
32 #else
33  #define EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN
34 #endif
35 
36 namespace Eigen {
37 
38 namespace internal {
39 
40 struct constructor_without_unaligned_array_assert {};
41 
46 template <typename T, int Size, int MatrixOrArrayOptions,
47  int Alignment = (MatrixOrArrayOptions&DontAlign) ? 0
48  : (((Size*sizeof(T))%16)==0) ? 16
49  : 0 >
50 struct plain_array
51 {
52  T array[Size];
53  plain_array() {}
54  plain_array(constructor_without_unaligned_array_assert) {}
55 };
56 
57 #ifdef EIGEN_DISABLE_UNALIGNED_ARRAY_ASSERT
58  #define EIGEN_MAKE_UNALIGNED_ARRAY_ASSERT(sizemask)
59 #else
60  #define EIGEN_MAKE_UNALIGNED_ARRAY_ASSERT(sizemask) \
61  eigen_assert((reinterpret_cast<size_t>(array) & sizemask) == 0 \
62  && "this assertion is explained here: " \
63  "http://eigen.tuxfamily.org/dox-devel/TopicUnalignedArrayAssert.html" \
64  " **** READ THIS WEB PAGE !!! ****");
65 #endif
66 
67 template <typename T, int Size, int MatrixOrArrayOptions>
68 struct plain_array<T, Size, MatrixOrArrayOptions, 16>
69 {
70  EIGEN_USER_ALIGN16 T array[Size];
71  plain_array() { EIGEN_MAKE_UNALIGNED_ARRAY_ASSERT(0xf) }
72  plain_array(constructor_without_unaligned_array_assert) {}
73 };
74 
75 template <typename T, int MatrixOrArrayOptions, int Alignment>
76 struct plain_array<T, 0, MatrixOrArrayOptions, Alignment>
77 {
78  EIGEN_USER_ALIGN16 T array[1];
79  plain_array() {}
80  plain_array(constructor_without_unaligned_array_assert) {}
81 };
82 
83 } // end namespace internal
84 
97 template<typename T, int Size, int _Rows, int _Cols, int _Options> class DenseStorage;
98 
99 // purely fixed-size matrix
100 template<typename T, int Size, int _Rows, int _Cols, int _Options> class DenseStorage
101 {
102  internal::plain_array<T,Size,_Options> m_data;
103  public:
104  inline explicit DenseStorage() {}
105  inline DenseStorage(internal::constructor_without_unaligned_array_assert)
106  : m_data(internal::constructor_without_unaligned_array_assert()) {}
108  inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); }
109  static inline DenseIndex rows(void) {return _Rows;}
110  static inline DenseIndex cols(void) {return _Cols;}
113  inline const T *data() const { return m_data.array; }
114  inline T *data() { return m_data.array; }
115 };
116 
117 // null matrix
118 template<typename T, int _Rows, int _Cols, int _Options> class DenseStorage<T, 0, _Rows, _Cols, _Options>
119 {
120  public:
121  inline explicit DenseStorage() {}
122  inline DenseStorage(internal::constructor_without_unaligned_array_assert) {}
123  inline DenseStorage(DenseIndex,DenseIndex,DenseIndex) {}
124  inline void swap(DenseStorage& ) {}
125  static inline DenseIndex rows(void) {return _Rows;}
126  static inline DenseIndex cols(void) {return _Cols;}
127  inline void conservativeResize(DenseIndex,DenseIndex,DenseIndex) {}
128  inline void resize(DenseIndex,DenseIndex,DenseIndex) {}
129  inline const T *data() const { return 0; }
130  inline T *data() { return 0; }
131 };
132 
133 // more specializations for null matrices; these are necessary to resolve ambiguities
134 template<typename T, int _Options> class DenseStorage<T, 0, Dynamic, Dynamic, _Options>
135 : public DenseStorage<T, 0, 0, 0, _Options> { };
136 
137 template<typename T, int _Rows, int _Options> class DenseStorage<T, 0, _Rows, Dynamic, _Options>
138 : public DenseStorage<T, 0, 0, 0, _Options> { };
139 
140 template<typename T, int _Cols, int _Options> class DenseStorage<T, 0, Dynamic, _Cols, _Options>
141 : public DenseStorage<T, 0, 0, 0, _Options> { };
142 
143 // dynamic-size matrix with fixed-size storage
144 template<typename T, int Size, int _Options> class DenseStorage<T, Size, Dynamic, Dynamic, _Options>
145 {
146  internal::plain_array<T,Size,_Options> m_data;
147  DenseIndex m_rows;
148  DenseIndex m_cols;
149  public:
150  inline explicit DenseStorage() : m_rows(0), m_cols(0) {}
151  inline DenseStorage(internal::constructor_without_unaligned_array_assert)
152  : m_data(internal::constructor_without_unaligned_array_assert()), m_rows(0), m_cols(0) {}
153  inline DenseStorage(DenseIndex, DenseIndex rows, DenseIndex cols) : m_rows(rows), m_cols(cols) {}
154  inline void swap(DenseStorage& other)
155  { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); std::swap(m_cols,other.m_cols); }
156  inline DenseIndex rows(void) const {return m_rows;}
157  inline DenseIndex cols(void) const {return m_cols;}
158  inline void conservativeResize(DenseIndex, DenseIndex rows, DenseIndex cols) { m_rows = rows; m_cols = cols; }
159  inline void resize(DenseIndex, DenseIndex rows, DenseIndex cols) { m_rows = rows; m_cols = cols; }
160  inline const T *data() const { return m_data.array; }
161  inline T *data() { return m_data.array; }
162 };
163 
164 // dynamic-size matrix with fixed-size storage and fixed width
165 template<typename T, int Size, int _Cols, int _Options> class DenseStorage<T, Size, Dynamic, _Cols, _Options>
166 {
167  internal::plain_array<T,Size,_Options> m_data;
168  DenseIndex m_rows;
169  public:
170  inline explicit DenseStorage() : m_rows(0) {}
171  inline DenseStorage(internal::constructor_without_unaligned_array_assert)
172  : m_data(internal::constructor_without_unaligned_array_assert()), m_rows(0) {}
173  inline DenseStorage(DenseIndex, DenseIndex rows, DenseIndex) : m_rows(rows) {}
174  inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); }
175  inline DenseIndex rows(void) const {return m_rows;}
176  inline DenseIndex cols(void) const {return _Cols;}
177  inline void conservativeResize(DenseIndex, DenseIndex rows, DenseIndex) { m_rows = rows; }
178  inline void resize(DenseIndex, DenseIndex rows, DenseIndex) { m_rows = rows; }
179  inline const T *data() const { return m_data.array; }
180  inline T *data() { return m_data.array; }
181 };
182 
183 // dynamic-size matrix with fixed-size storage and fixed height
184 template<typename T, int Size, int _Rows, int _Options> class DenseStorage<T, Size, _Rows, Dynamic, _Options>
185 {
186  internal::plain_array<T,Size,_Options> m_data;
187  DenseIndex m_cols;
188  public:
189  inline explicit DenseStorage() : m_cols(0) {}
190  inline DenseStorage(internal::constructor_without_unaligned_array_assert)
191  : m_data(internal::constructor_without_unaligned_array_assert()), m_cols(0) {}
192  inline DenseStorage(DenseIndex, DenseIndex, DenseIndex cols) : m_cols(cols) {}
193  inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_cols,other.m_cols); }
194  inline DenseIndex rows(void) const {return _Rows;}
195  inline DenseIndex cols(void) const {return m_cols;}
196  inline void conservativeResize(DenseIndex, DenseIndex, DenseIndex cols) { m_cols = cols; }
197  inline void resize(DenseIndex, DenseIndex, DenseIndex cols) { m_cols = cols; }
198  inline const T *data() const { return m_data.array; }
199  inline T *data() { return m_data.array; }
200 };
201 
202 // purely dynamic matrix.
203 template<typename T, int _Options> class DenseStorage<T, Dynamic, Dynamic, Dynamic, _Options>
204 {
205  T *m_data;
206  DenseIndex m_rows;
207  DenseIndex m_cols;
208  public:
209  inline explicit DenseStorage() : m_data(0), m_rows(0), m_cols(0) {}
210  inline DenseStorage(internal::constructor_without_unaligned_array_assert)
211  : m_data(0), m_rows(0), m_cols(0) {}
212  inline DenseStorage(DenseIndex size, DenseIndex rows, DenseIndex cols)
213  : m_data(internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size)), m_rows(rows), m_cols(cols)
215  inline ~DenseStorage() { internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, m_rows*m_cols); }
216  inline void swap(DenseStorage& other)
217  { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); std::swap(m_cols,other.m_cols); }
218  inline DenseIndex rows(void) const {return m_rows;}
219  inline DenseIndex cols(void) const {return m_cols;}
220  inline void conservativeResize(DenseIndex size, DenseIndex rows, DenseIndex cols)
221  {
222  m_data = internal::conditional_aligned_realloc_new_auto<T,(_Options&DontAlign)==0>(m_data, size, m_rows*m_cols);
223  m_rows = rows;
224  m_cols = cols;
225  }
226  void resize(DenseIndex size, DenseIndex rows, DenseIndex cols)
227  {
228  if(size != m_rows*m_cols)
229  {
230  internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, m_rows*m_cols);
231  if (size)
232  m_data = internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size);
233  else
234  m_data = 0;
236  }
237  m_rows = rows;
238  m_cols = cols;
239  }
240  inline const T *data() const { return m_data; }
241  inline T *data() { return m_data; }
242 };
243 
244 // matrix with dynamic width and fixed height (so that matrix has dynamic size).
245 template<typename T, int _Rows, int _Options> class DenseStorage<T, Dynamic, _Rows, Dynamic, _Options>
246 {
247  T *m_data;
248  DenseIndex m_cols;
249  public:
250  inline explicit DenseStorage() : m_data(0), m_cols(0) {}
251  inline DenseStorage(internal::constructor_without_unaligned_array_assert) : m_data(0), m_cols(0) {}
252  inline DenseStorage(DenseIndex size, DenseIndex, DenseIndex cols) : m_data(internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size)), m_cols(cols)
254  inline ~DenseStorage() { internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Rows*m_cols); }
255  inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_cols,other.m_cols); }
256  static inline DenseIndex rows(void) {return _Rows;}
257  inline DenseIndex cols(void) const {return m_cols;}
258  inline void conservativeResize(DenseIndex size, DenseIndex, DenseIndex cols)
259  {
260  m_data = internal::conditional_aligned_realloc_new_auto<T,(_Options&DontAlign)==0>(m_data, size, _Rows*m_cols);
261  m_cols = cols;
262  }
263  EIGEN_STRONG_INLINE void resize(DenseIndex size, DenseIndex, DenseIndex cols)
264  {
265  if(size != _Rows*m_cols)
266  {
267  internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Rows*m_cols);
268  if (size)
269  m_data = internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size);
270  else
271  m_data = 0;
273  }
274  m_cols = cols;
275  }
276  inline const T *data() const { return m_data; }
277  inline T *data() { return m_data; }
278 };
279 
280 // matrix with dynamic height and fixed width (so that matrix has dynamic size).
281 template<typename T, int _Cols, int _Options> class DenseStorage<T, Dynamic, Dynamic, _Cols, _Options>
282 {
283  T *m_data;
284  DenseIndex m_rows;
285  public:
286  inline explicit DenseStorage() : m_data(0), m_rows(0) {}
287  inline DenseStorage(internal::constructor_without_unaligned_array_assert) : m_data(0), m_rows(0) {}
288  inline DenseStorage(DenseIndex size, DenseIndex rows, DenseIndex) : m_data(internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size)), m_rows(rows)
290  inline ~DenseStorage() { internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Cols*m_rows); }
291  inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); }
292  inline DenseIndex rows(void) const {return m_rows;}
293  static inline DenseIndex cols(void) {return _Cols;}
294  inline void conservativeResize(DenseIndex size, DenseIndex rows, DenseIndex)
295  {
296  m_data = internal::conditional_aligned_realloc_new_auto<T,(_Options&DontAlign)==0>(m_data, size, m_rows*_Cols);
297  m_rows = rows;
298  }
299  EIGEN_STRONG_INLINE void resize(DenseIndex size, DenseIndex rows, DenseIndex)
300  {
301  if(size != m_rows*_Cols)
302  {
303  internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Cols*m_rows);
304  if (size)
305  m_data = internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size);
306  else
307  m_data = 0;
309  }
310  m_rows = rows;
311  }
312  inline const T *data() const { return m_data; }
313  inline T *data() { return m_data; }
314 };
315 
316 } // end namespace Eigen
317 
318 #endif // EIGEN_MATRIX_H