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// This Source Code Form is subject to the terms of the Mozilla
9// Public License v. 2.0. If a copy of the MPL was not distributed
10// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
11
12#ifndef EIGEN_MATRIXSTORAGE_H
13#define EIGEN_MATRIXSTORAGE_H
14
15#ifdef EIGEN_DENSE_STORAGE_CTOR_PLUGIN
16  #define EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN EIGEN_DENSE_STORAGE_CTOR_PLUGIN;
17#else
18  #define EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN
19#endif
20
21namespace Eigen {
22
23namespace internal {
24
25struct constructor_without_unaligned_array_assert {};
26
27template<typename T, int Size> void check_static_allocation_size()
28{
29  // if EIGEN_STACK_ALLOCATION_LIMIT is defined to 0, then no limit
30  #if EIGEN_STACK_ALLOCATION_LIMIT
31  EIGEN_STATIC_ASSERT(Size * sizeof(T) <= EIGEN_STACK_ALLOCATION_LIMIT, OBJECT_ALLOCATED_ON_STACK_IS_TOO_BIG);
32  #endif
33}
34
35/** \internal
36  * Static array. If the MatrixOrArrayOptions require auto-alignment, the array will be automatically aligned:
37  * to 16 bytes boundary if the total size is a multiple of 16 bytes.
38  */
39template <typename T, int Size, int MatrixOrArrayOptions,
40          int Alignment = (MatrixOrArrayOptions&DontAlign) ? 0
41                        : (((Size*sizeof(T))%16)==0) ? 16
42                        : 0 >
43struct plain_array
44{
45  T array[Size];
46
47  plain_array()
48  {
49    check_static_allocation_size<T,Size>();
50  }
51
52  plain_array(constructor_without_unaligned_array_assert)
53  {
54    check_static_allocation_size<T,Size>();
55  }
56};
57
58#if defined(EIGEN_DISABLE_UNALIGNED_ARRAY_ASSERT)
59  #define EIGEN_MAKE_UNALIGNED_ARRAY_ASSERT(sizemask)
60#elif EIGEN_GNUC_AT_LEAST(4,7)
61  // GCC 4.7 is too aggressive in its optimizations and remove the alignement test based on the fact the array is declared to be aligned.
62  // See this bug report: http://gcc.gnu.org/bugzilla/show_bug.cgi?id=53900
63  // Hiding the origin of the array pointer behind a function argument seems to do the trick even if the function is inlined:
64  template<typename PtrType>
65  EIGEN_ALWAYS_INLINE PtrType eigen_unaligned_array_assert_workaround_gcc47(PtrType array) { return array; }
66  #define EIGEN_MAKE_UNALIGNED_ARRAY_ASSERT(sizemask) \
67    eigen_assert((reinterpret_cast<size_t>(eigen_unaligned_array_assert_workaround_gcc47(array)) & sizemask) == 0 \
68              && "this assertion is explained here: " \
69              "http://eigen.tuxfamily.org/dox-devel/group__TopicUnalignedArrayAssert.html" \
70              " **** READ THIS WEB PAGE !!! ****");
71#else
72  #define EIGEN_MAKE_UNALIGNED_ARRAY_ASSERT(sizemask) \
73    eigen_assert((reinterpret_cast<size_t>(array) & sizemask) == 0 \
74              && "this assertion is explained here: " \
75              "http://eigen.tuxfamily.org/dox-devel/group__TopicUnalignedArrayAssert.html" \
76              " **** READ THIS WEB PAGE !!! ****");
77#endif
78
79template <typename T, int Size, int MatrixOrArrayOptions>
80struct plain_array<T, Size, MatrixOrArrayOptions, 16>
81{
82  EIGEN_USER_ALIGN16 T array[Size];
83
84  plain_array()
85  {
86    EIGEN_MAKE_UNALIGNED_ARRAY_ASSERT(0xf);
87    check_static_allocation_size<T,Size>();
88  }
89
90  plain_array(constructor_without_unaligned_array_assert)
91  {
92    check_static_allocation_size<T,Size>();
93  }
94};
95
96template <typename T, int MatrixOrArrayOptions, int Alignment>
97struct plain_array<T, 0, MatrixOrArrayOptions, Alignment>
98{
99  EIGEN_USER_ALIGN16 T array[1];
100  plain_array() {}
101  plain_array(constructor_without_unaligned_array_assert) {}
102};
103
104} // end namespace internal
105
106/** \internal
107  *
108  * \class DenseStorage
109  * \ingroup Core_Module
110  *
111  * \brief Stores the data of a matrix
112  *
113  * This class stores the data of fixed-size, dynamic-size or mixed matrices
114  * in a way as compact as possible.
115  *
116  * \sa Matrix
117  */
118template<typename T, int Size, int _Rows, int _Cols, int _Options> class DenseStorage;
119
120// purely fixed-size matrix
121template<typename T, int Size, int _Rows, int _Cols, int _Options> class DenseStorage
122{
123    internal::plain_array<T,Size,_Options> m_data;
124  public:
125    inline DenseStorage() {}
126    inline DenseStorage(internal::constructor_without_unaligned_array_assert)
127      : m_data(internal::constructor_without_unaligned_array_assert()) {}
128    inline DenseStorage(DenseIndex,DenseIndex,DenseIndex) {}
129    inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); }
130    static inline DenseIndex rows(void) {return _Rows;}
131    static inline DenseIndex cols(void) {return _Cols;}
132    inline void conservativeResize(DenseIndex,DenseIndex,DenseIndex) {}
133    inline void resize(DenseIndex,DenseIndex,DenseIndex) {}
134    inline const T *data() const { return m_data.array; }
135    inline T *data() { return m_data.array; }
136};
137
138// null matrix
139template<typename T, int _Rows, int _Cols, int _Options> class DenseStorage<T, 0, _Rows, _Cols, _Options>
140{
141  public:
142    inline DenseStorage() {}
143    inline DenseStorage(internal::constructor_without_unaligned_array_assert) {}
144    inline DenseStorage(DenseIndex,DenseIndex,DenseIndex) {}
145    inline void swap(DenseStorage& ) {}
146    static inline DenseIndex rows(void) {return _Rows;}
147    static inline DenseIndex cols(void) {return _Cols;}
148    inline void conservativeResize(DenseIndex,DenseIndex,DenseIndex) {}
149    inline void resize(DenseIndex,DenseIndex,DenseIndex) {}
150    inline const T *data() const { return 0; }
151    inline T *data() { return 0; }
152};
153
154// more specializations for null matrices; these are necessary to resolve ambiguities
155template<typename T, int _Options> class DenseStorage<T, 0, Dynamic, Dynamic, _Options>
156: public DenseStorage<T, 0, 0, 0, _Options> { };
157
158template<typename T, int _Rows, int _Options> class DenseStorage<T, 0, _Rows, Dynamic, _Options>
159: public DenseStorage<T, 0, 0, 0, _Options> { };
160
161template<typename T, int _Cols, int _Options> class DenseStorage<T, 0, Dynamic, _Cols, _Options>
162: public DenseStorage<T, 0, 0, 0, _Options> { };
163
164// dynamic-size matrix with fixed-size storage
165template<typename T, int Size, int _Options> class DenseStorage<T, Size, Dynamic, Dynamic, _Options>
166{
167    internal::plain_array<T,Size,_Options> m_data;
168    DenseIndex m_rows;
169    DenseIndex m_cols;
170  public:
171    inline DenseStorage() : m_rows(0), m_cols(0) {}
172    inline DenseStorage(internal::constructor_without_unaligned_array_assert)
173      : m_data(internal::constructor_without_unaligned_array_assert()), m_rows(0), m_cols(0) {}
174    inline DenseStorage(DenseIndex, DenseIndex nbRows, DenseIndex nbCols) : m_rows(nbRows), m_cols(nbCols) {}
175    inline void swap(DenseStorage& other)
176    { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); std::swap(m_cols,other.m_cols); }
177    inline DenseIndex rows() const {return m_rows;}
178    inline DenseIndex cols() const {return m_cols;}
179    inline void conservativeResize(DenseIndex, DenseIndex nbRows, DenseIndex nbCols) { m_rows = nbRows; m_cols = nbCols; }
180    inline void resize(DenseIndex, DenseIndex nbRows, DenseIndex nbCols) { m_rows = nbRows; m_cols = nbCols; }
181    inline const T *data() const { return m_data.array; }
182    inline T *data() { return m_data.array; }
183};
184
185// dynamic-size matrix with fixed-size storage and fixed width
186template<typename T, int Size, int _Cols, int _Options> class DenseStorage<T, Size, Dynamic, _Cols, _Options>
187{
188    internal::plain_array<T,Size,_Options> m_data;
189    DenseIndex m_rows;
190  public:
191    inline DenseStorage() : m_rows(0) {}
192    inline DenseStorage(internal::constructor_without_unaligned_array_assert)
193      : m_data(internal::constructor_without_unaligned_array_assert()), m_rows(0) {}
194    inline DenseStorage(DenseIndex, DenseIndex nbRows, DenseIndex) : m_rows(nbRows) {}
195    inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); }
196    inline DenseIndex rows(void) const {return m_rows;}
197    inline DenseIndex cols(void) const {return _Cols;}
198    inline void conservativeResize(DenseIndex, DenseIndex nbRows, DenseIndex) { m_rows = nbRows; }
199    inline void resize(DenseIndex, DenseIndex nbRows, DenseIndex) { m_rows = nbRows; }
200    inline const T *data() const { return m_data.array; }
201    inline T *data() { return m_data.array; }
202};
203
204// dynamic-size matrix with fixed-size storage and fixed height
205template<typename T, int Size, int _Rows, int _Options> class DenseStorage<T, Size, _Rows, Dynamic, _Options>
206{
207    internal::plain_array<T,Size,_Options> m_data;
208    DenseIndex m_cols;
209  public:
210    inline DenseStorage() : m_cols(0) {}
211    inline DenseStorage(internal::constructor_without_unaligned_array_assert)
212      : m_data(internal::constructor_without_unaligned_array_assert()), m_cols(0) {}
213    inline DenseStorage(DenseIndex, DenseIndex, DenseIndex nbCols) : m_cols(nbCols) {}
214    inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_cols,other.m_cols); }
215    inline DenseIndex rows(void) const {return _Rows;}
216    inline DenseIndex cols(void) const {return m_cols;}
217    inline void conservativeResize(DenseIndex, DenseIndex, DenseIndex nbCols) { m_cols = nbCols; }
218    inline void resize(DenseIndex, DenseIndex, DenseIndex nbCols) { m_cols = nbCols; }
219    inline const T *data() const { return m_data.array; }
220    inline T *data() { return m_data.array; }
221};
222
223// purely dynamic matrix.
224template<typename T, int _Options> class DenseStorage<T, Dynamic, Dynamic, Dynamic, _Options>
225{
226    T *m_data;
227    DenseIndex m_rows;
228    DenseIndex m_cols;
229  public:
230    inline DenseStorage() : m_data(0), m_rows(0), m_cols(0) {}
231    inline DenseStorage(internal::constructor_without_unaligned_array_assert)
232       : m_data(0), m_rows(0), m_cols(0) {}
233    inline DenseStorage(DenseIndex size, DenseIndex nbRows, DenseIndex nbCols)
234      : m_data(internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size)), m_rows(nbRows), m_cols(nbCols)
235    { EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN }
236    inline ~DenseStorage() { internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, m_rows*m_cols); }
237    inline void swap(DenseStorage& other)
238    { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); std::swap(m_cols,other.m_cols); }
239    inline DenseIndex rows(void) const {return m_rows;}
240    inline DenseIndex cols(void) const {return m_cols;}
241    inline void conservativeResize(DenseIndex size, DenseIndex nbRows, DenseIndex nbCols)
242    {
243      m_data = internal::conditional_aligned_realloc_new_auto<T,(_Options&DontAlign)==0>(m_data, size, m_rows*m_cols);
244      m_rows = nbRows;
245      m_cols = nbCols;
246    }
247    void resize(DenseIndex size, DenseIndex nbRows, DenseIndex nbCols)
248    {
249      if(size != m_rows*m_cols)
250      {
251        internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, m_rows*m_cols);
252        if (size)
253          m_data = internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size);
254        else
255          m_data = 0;
256        EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN
257      }
258      m_rows = nbRows;
259      m_cols = nbCols;
260    }
261    inline const T *data() const { return m_data; }
262    inline T *data() { return m_data; }
263};
264
265// matrix with dynamic width and fixed height (so that matrix has dynamic size).
266template<typename T, int _Rows, int _Options> class DenseStorage<T, Dynamic, _Rows, Dynamic, _Options>
267{
268    T *m_data;
269    DenseIndex m_cols;
270  public:
271    inline DenseStorage() : m_data(0), m_cols(0) {}
272    inline DenseStorage(internal::constructor_without_unaligned_array_assert) : m_data(0), m_cols(0) {}
273    inline DenseStorage(DenseIndex size, DenseIndex, DenseIndex nbCols) : m_data(internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size)), m_cols(nbCols)
274    { EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN }
275    inline ~DenseStorage() { internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Rows*m_cols); }
276    inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_cols,other.m_cols); }
277    static inline DenseIndex rows(void) {return _Rows;}
278    inline DenseIndex cols(void) const {return m_cols;}
279    inline void conservativeResize(DenseIndex size, DenseIndex, DenseIndex nbCols)
280    {
281      m_data = internal::conditional_aligned_realloc_new_auto<T,(_Options&DontAlign)==0>(m_data, size, _Rows*m_cols);
282      m_cols = nbCols;
283    }
284    EIGEN_STRONG_INLINE void resize(DenseIndex size, DenseIndex, DenseIndex nbCols)
285    {
286      if(size != _Rows*m_cols)
287      {
288        internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Rows*m_cols);
289        if (size)
290          m_data = internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size);
291        else
292          m_data = 0;
293        EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN
294      }
295      m_cols = nbCols;
296    }
297    inline const T *data() const { return m_data; }
298    inline T *data() { return m_data; }
299};
300
301// matrix with dynamic height and fixed width (so that matrix has dynamic size).
302template<typename T, int _Cols, int _Options> class DenseStorage<T, Dynamic, Dynamic, _Cols, _Options>
303{
304    T *m_data;
305    DenseIndex m_rows;
306  public:
307    inline DenseStorage() : m_data(0), m_rows(0) {}
308    inline DenseStorage(internal::constructor_without_unaligned_array_assert) : m_data(0), m_rows(0) {}
309    inline DenseStorage(DenseIndex size, DenseIndex nbRows, DenseIndex) : m_data(internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size)), m_rows(nbRows)
310    { EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN }
311    inline ~DenseStorage() { internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Cols*m_rows); }
312    inline void swap(DenseStorage& other) { std::swap(m_data,other.m_data); std::swap(m_rows,other.m_rows); }
313    inline DenseIndex rows(void) const {return m_rows;}
314    static inline DenseIndex cols(void) {return _Cols;}
315    inline void conservativeResize(DenseIndex size, DenseIndex nbRows, DenseIndex)
316    {
317      m_data = internal::conditional_aligned_realloc_new_auto<T,(_Options&DontAlign)==0>(m_data, size, m_rows*_Cols);
318      m_rows = nbRows;
319    }
320    EIGEN_STRONG_INLINE void resize(DenseIndex size, DenseIndex nbRows, DenseIndex)
321    {
322      if(size != m_rows*_Cols)
323      {
324        internal::conditional_aligned_delete_auto<T,(_Options&DontAlign)==0>(m_data, _Cols*m_rows);
325        if (size)
326          m_data = internal::conditional_aligned_new_auto<T,(_Options&DontAlign)==0>(size);
327        else
328          m_data = 0;
329        EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN
330      }
331      m_rows = nbRows;
332    }
333    inline const T *data() const { return m_data; }
334    inline T *data() { return m_data; }
335};
336
337} // end namespace Eigen
338
339#endif // EIGEN_MATRIX_H
340