1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2009-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
5//
6// This Source Code Form is subject to the terms of the Mozilla
7// Public License v. 2.0. If a copy of the MPL was not distributed
8// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9
10#ifndef EIGEN_ARRAYWRAPPER_H
11#define EIGEN_ARRAYWRAPPER_H
12
13namespace Eigen {
14
15/** \class ArrayWrapper
16  * \ingroup Core_Module
17  *
18  * \brief Expression of a mathematical vector or matrix as an array object
19  *
20  * This class is the return type of MatrixBase::array(), and most of the time
21  * this is the only way it is use.
22  *
23  * \sa MatrixBase::array(), class MatrixWrapper
24  */
25
26namespace internal {
27template<typename ExpressionType>
28struct traits<ArrayWrapper<ExpressionType> >
29  : public traits<typename remove_all<typename ExpressionType::Nested>::type >
30{
31  typedef ArrayXpr XprKind;
32};
33}
34
35template<typename ExpressionType>
36class ArrayWrapper : public ArrayBase<ArrayWrapper<ExpressionType> >
37{
38  public:
39    typedef ArrayBase<ArrayWrapper> Base;
40    EIGEN_DENSE_PUBLIC_INTERFACE(ArrayWrapper)
41    EIGEN_INHERIT_ASSIGNMENT_OPERATORS(ArrayWrapper)
42
43    typedef typename internal::conditional<
44                       internal::is_lvalue<ExpressionType>::value,
45                       Scalar,
46                       const Scalar
47                     >::type ScalarWithConstIfNotLvalue;
48
49    typedef typename internal::nested<ExpressionType>::type NestedExpressionType;
50
51    inline ArrayWrapper(ExpressionType& matrix) : m_expression(matrix) {}
52
53    inline Index rows() const { return m_expression.rows(); }
54    inline Index cols() const { return m_expression.cols(); }
55    inline Index outerStride() const { return m_expression.outerStride(); }
56    inline Index innerStride() const { return m_expression.innerStride(); }
57
58    inline ScalarWithConstIfNotLvalue* data() { return m_expression.const_cast_derived().data(); }
59    inline const Scalar* data() const { return m_expression.data(); }
60
61    inline CoeffReturnType coeff(Index rowId, Index colId) const
62    {
63      return m_expression.coeff(rowId, colId);
64    }
65
66    inline Scalar& coeffRef(Index rowId, Index colId)
67    {
68      return m_expression.const_cast_derived().coeffRef(rowId, colId);
69    }
70
71    inline const Scalar& coeffRef(Index rowId, Index colId) const
72    {
73      return m_expression.const_cast_derived().coeffRef(rowId, colId);
74    }
75
76    inline CoeffReturnType coeff(Index index) const
77    {
78      return m_expression.coeff(index);
79    }
80
81    inline Scalar& coeffRef(Index index)
82    {
83      return m_expression.const_cast_derived().coeffRef(index);
84    }
85
86    inline const Scalar& coeffRef(Index index) const
87    {
88      return m_expression.const_cast_derived().coeffRef(index);
89    }
90
91    template<int LoadMode>
92    inline const PacketScalar packet(Index rowId, Index colId) const
93    {
94      return m_expression.template packet<LoadMode>(rowId, colId);
95    }
96
97    template<int LoadMode>
98    inline void writePacket(Index rowId, Index colId, const PacketScalar& val)
99    {
100      m_expression.const_cast_derived().template writePacket<LoadMode>(rowId, colId, val);
101    }
102
103    template<int LoadMode>
104    inline const PacketScalar packet(Index index) const
105    {
106      return m_expression.template packet<LoadMode>(index);
107    }
108
109    template<int LoadMode>
110    inline void writePacket(Index index, const PacketScalar& val)
111    {
112      m_expression.const_cast_derived().template writePacket<LoadMode>(index, val);
113    }
114
115    template<typename Dest>
116    inline void evalTo(Dest& dst) const { dst = m_expression; }
117
118    const typename internal::remove_all<NestedExpressionType>::type&
119    nestedExpression() const
120    {
121      return m_expression;
122    }
123
124    /** Forwards the resizing request to the nested expression
125      * \sa DenseBase::resize(Index)  */
126    void resize(Index newSize) { m_expression.const_cast_derived().resize(newSize); }
127    /** Forwards the resizing request to the nested expression
128      * \sa DenseBase::resize(Index,Index)*/
129    void resize(Index nbRows, Index nbCols) { m_expression.const_cast_derived().resize(nbRows,nbCols); }
130
131  protected:
132    NestedExpressionType m_expression;
133};
134
135/** \class MatrixWrapper
136  * \ingroup Core_Module
137  *
138  * \brief Expression of an array as a mathematical vector or matrix
139  *
140  * This class is the return type of ArrayBase::matrix(), and most of the time
141  * this is the only way it is use.
142  *
143  * \sa MatrixBase::matrix(), class ArrayWrapper
144  */
145
146namespace internal {
147template<typename ExpressionType>
148struct traits<MatrixWrapper<ExpressionType> >
149 : public traits<typename remove_all<typename ExpressionType::Nested>::type >
150{
151  typedef MatrixXpr XprKind;
152};
153}
154
155template<typename ExpressionType>
156class MatrixWrapper : public MatrixBase<MatrixWrapper<ExpressionType> >
157{
158  public:
159    typedef MatrixBase<MatrixWrapper<ExpressionType> > Base;
160    EIGEN_DENSE_PUBLIC_INTERFACE(MatrixWrapper)
161    EIGEN_INHERIT_ASSIGNMENT_OPERATORS(MatrixWrapper)
162
163    typedef typename internal::conditional<
164                       internal::is_lvalue<ExpressionType>::value,
165                       Scalar,
166                       const Scalar
167                     >::type ScalarWithConstIfNotLvalue;
168
169    typedef typename internal::nested<ExpressionType>::type NestedExpressionType;
170
171    inline MatrixWrapper(ExpressionType& a_matrix) : m_expression(a_matrix) {}
172
173    inline Index rows() const { return m_expression.rows(); }
174    inline Index cols() const { return m_expression.cols(); }
175    inline Index outerStride() const { return m_expression.outerStride(); }
176    inline Index innerStride() const { return m_expression.innerStride(); }
177
178    inline ScalarWithConstIfNotLvalue* data() { return m_expression.const_cast_derived().data(); }
179    inline const Scalar* data() const { return m_expression.data(); }
180
181    inline CoeffReturnType coeff(Index rowId, Index colId) const
182    {
183      return m_expression.coeff(rowId, colId);
184    }
185
186    inline Scalar& coeffRef(Index rowId, Index colId)
187    {
188      return m_expression.const_cast_derived().coeffRef(rowId, colId);
189    }
190
191    inline const Scalar& coeffRef(Index rowId, Index colId) const
192    {
193      return m_expression.derived().coeffRef(rowId, colId);
194    }
195
196    inline CoeffReturnType coeff(Index index) const
197    {
198      return m_expression.coeff(index);
199    }
200
201    inline Scalar& coeffRef(Index index)
202    {
203      return m_expression.const_cast_derived().coeffRef(index);
204    }
205
206    inline const Scalar& coeffRef(Index index) const
207    {
208      return m_expression.const_cast_derived().coeffRef(index);
209    }
210
211    template<int LoadMode>
212    inline const PacketScalar packet(Index rowId, Index colId) const
213    {
214      return m_expression.template packet<LoadMode>(rowId, colId);
215    }
216
217    template<int LoadMode>
218    inline void writePacket(Index rowId, Index colId, const PacketScalar& val)
219    {
220      m_expression.const_cast_derived().template writePacket<LoadMode>(rowId, colId, val);
221    }
222
223    template<int LoadMode>
224    inline const PacketScalar packet(Index index) const
225    {
226      return m_expression.template packet<LoadMode>(index);
227    }
228
229    template<int LoadMode>
230    inline void writePacket(Index index, const PacketScalar& val)
231    {
232      m_expression.const_cast_derived().template writePacket<LoadMode>(index, val);
233    }
234
235    const typename internal::remove_all<NestedExpressionType>::type&
236    nestedExpression() const
237    {
238      return m_expression;
239    }
240
241    /** Forwards the resizing request to the nested expression
242      * \sa DenseBase::resize(Index)  */
243    void resize(Index newSize) { m_expression.const_cast_derived().resize(newSize); }
244    /** Forwards the resizing request to the nested expression
245      * \sa DenseBase::resize(Index,Index)*/
246    void resize(Index nbRows, Index nbCols) { m_expression.const_cast_derived().resize(nbRows,nbCols); }
247
248  protected:
249    NestedExpressionType m_expression;
250};
251
252} // end namespace Eigen
253
254#endif // EIGEN_ARRAYWRAPPER_H
255