Eigen  3.3.90 (mercurial changeset 493691b29be1)
IndexedViewMethods.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2017 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 #if !defined(EIGEN_PARSED_BY_DOXYGEN)
11 
12 // This file is automatically included twice to generate const and non-const versions
13 
14 #ifndef EIGEN_INDEXED_VIEW_METHOD_2ND_PASS
15 #define EIGEN_INDEXED_VIEW_METHOD_CONST const
16 #define EIGEN_INDEXED_VIEW_METHOD_TYPE ConstIndexedViewType
17 #else
18 #define EIGEN_INDEXED_VIEW_METHOD_CONST
19 #define EIGEN_INDEXED_VIEW_METHOD_TYPE IndexedViewType
20 #endif
21 
22 #ifndef EIGEN_INDEXED_VIEW_METHOD_2ND_PASS
23 protected:
24 
25 // define some aliases to ease readability
26 
27 template<typename Indices>
28 struct IvcRowType : public internal::IndexedViewCompatibleType<Indices,RowsAtCompileTime> {};
29 
30 template<typename Indices>
31 struct IvcColType : public internal::IndexedViewCompatibleType<Indices,ColsAtCompileTime> {};
32 
33 template<typename Indices>
34 struct IvcType : public internal::IndexedViewCompatibleType<Indices,SizeAtCompileTime> {};
35 
36 typedef typename internal::IndexedViewCompatibleType<Index,1>::type IvcIndex;
37 
38 template<typename Indices>
39 typename IvcRowType<Indices>::type
40 ivcRow(const Indices& indices) const {
41  return internal::makeIndexedViewCompatible(indices, internal::variable_if_dynamic<Index,RowsAtCompileTime>(derived().rows()),Specialized);
42 }
43 
44 template<typename Indices>
45 typename IvcColType<Indices>::type
46 ivcCol(const Indices& indices) const {
47  return internal::makeIndexedViewCompatible(indices, internal::variable_if_dynamic<Index,ColsAtCompileTime>(derived().cols()),Specialized);
48 }
49 
50 template<typename Indices>
51 typename IvcColType<Indices>::type
52 ivcSize(const Indices& indices) const {
53  return internal::makeIndexedViewCompatible(indices, internal::variable_if_dynamic<Index,SizeAtCompileTime>(derived().size()),Specialized);
54 }
55 
56 template<typename RowIndices, typename ColIndices>
57 struct valid_indexed_view_overload {
58  enum { value = !(internal::is_valid_index_type<RowIndices>::value && internal::is_valid_index_type<ColIndices>::value) };
59 };
60 
61 public:
62 
63 #endif
64 
65 template<typename RowIndices, typename ColIndices>
66 struct EIGEN_INDEXED_VIEW_METHOD_TYPE {
67  typedef IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,
68  typename IvcRowType<RowIndices>::type,
69  typename IvcColType<ColIndices>::type> type;
70 };
71 
72 // This is the generic version
73 
74 template<typename RowIndices, typename ColIndices>
75 typename internal::enable_if<valid_indexed_view_overload<RowIndices,ColIndices>::value
76  && internal::traits<typename EIGEN_INDEXED_VIEW_METHOD_TYPE<RowIndices,ColIndices>::type>::ReturnAsIndexedView,
77  typename EIGEN_INDEXED_VIEW_METHOD_TYPE<RowIndices,ColIndices>::type >::type
78 operator()(const RowIndices& rowIndices, const ColIndices& colIndices) EIGEN_INDEXED_VIEW_METHOD_CONST
79 {
80  return typename EIGEN_INDEXED_VIEW_METHOD_TYPE<RowIndices,ColIndices>::type
81  (derived(), ivcRow(rowIndices), ivcCol(colIndices));
82 }
83 
84 // The following overload returns a Block<> object
85 
86 template<typename RowIndices, typename ColIndices>
87 typename internal::enable_if<valid_indexed_view_overload<RowIndices,ColIndices>::value
88  && internal::traits<typename EIGEN_INDEXED_VIEW_METHOD_TYPE<RowIndices,ColIndices>::type>::ReturnAsBlock,
89  typename internal::traits<typename EIGEN_INDEXED_VIEW_METHOD_TYPE<RowIndices,ColIndices>::type>::BlockType>::type
90 operator()(const RowIndices& rowIndices, const ColIndices& colIndices) EIGEN_INDEXED_VIEW_METHOD_CONST
91 {
92  typedef typename internal::traits<typename EIGEN_INDEXED_VIEW_METHOD_TYPE<RowIndices,ColIndices>::type>::BlockType BlockType;
93  typename IvcRowType<RowIndices>::type actualRowIndices = ivcRow(rowIndices);
94  typename IvcColType<ColIndices>::type actualColIndices = ivcCol(colIndices);
95  return BlockType(derived(),
96  internal::first(actualRowIndices),
97  internal::first(actualColIndices),
98  internal::size(actualRowIndices),
99  internal::size(actualColIndices));
100 }
101 
102 // The following overload returns a Scalar
103 
104 template<typename RowIndices, typename ColIndices>
105 typename internal::enable_if<valid_indexed_view_overload<RowIndices,ColIndices>::value
106  && internal::traits<typename EIGEN_INDEXED_VIEW_METHOD_TYPE<RowIndices,ColIndices>::type>::ReturnAsScalar,
107  CoeffReturnType >::type
108 operator()(const RowIndices& rowIndices, const ColIndices& colIndices) EIGEN_INDEXED_VIEW_METHOD_CONST
109 {
110  return Base::operator()(internal::eval_expr_given_size(rowIndices,rows()),internal::eval_expr_given_size(colIndices,cols()));
111 }
112 
113 #if EIGEN_HAS_STATIC_ARRAY_TEMPLATE
114 
115 // The folowing three overloads are needed to handle raw Index[N] arrays.
116 
117 template<typename RowIndicesT, std::size_t RowIndicesN, typename ColIndices>
118 IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,const RowIndicesT (&)[RowIndicesN],typename IvcColType<ColIndices>::type>
119 operator()(const RowIndicesT (&rowIndices)[RowIndicesN], const ColIndices& colIndices) EIGEN_INDEXED_VIEW_METHOD_CONST
120 {
121  return IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,const RowIndicesT (&)[RowIndicesN],typename IvcColType<ColIndices>::type>
122  (derived(), rowIndices, ivcCol(colIndices));
123 }
124 
125 template<typename RowIndices, typename ColIndicesT, std::size_t ColIndicesN>
126 IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,typename IvcRowType<RowIndices>::type, const ColIndicesT (&)[ColIndicesN]>
127 operator()(const RowIndices& rowIndices, const ColIndicesT (&colIndices)[ColIndicesN]) EIGEN_INDEXED_VIEW_METHOD_CONST
128 {
129  return IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,typename IvcRowType<RowIndices>::type,const ColIndicesT (&)[ColIndicesN]>
130  (derived(), ivcRow(rowIndices), colIndices);
131 }
132 
133 template<typename RowIndicesT, std::size_t RowIndicesN, typename ColIndicesT, std::size_t ColIndicesN>
134 IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,const RowIndicesT (&)[RowIndicesN], const ColIndicesT (&)[ColIndicesN]>
135 operator()(const RowIndicesT (&rowIndices)[RowIndicesN], const ColIndicesT (&colIndices)[ColIndicesN]) EIGEN_INDEXED_VIEW_METHOD_CONST
136 {
137  return IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,const RowIndicesT (&)[RowIndicesN],const ColIndicesT (&)[ColIndicesN]>
138  (derived(), rowIndices, colIndices);
139 }
140 
141 #endif // EIGEN_HAS_STATIC_ARRAY_TEMPLATE
142 
143 // Overloads for 1D vectors/arrays
144 
145 template<typename Indices>
146 typename internal::enable_if<
147  IsRowMajor && (!(internal::get_compile_time_incr<typename IvcType<Indices>::type>::value==1 || internal::is_valid_index_type<Indices>::value)),
148  IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,IvcIndex,typename IvcType<Indices>::type> >::type
149 operator()(const Indices& indices) EIGEN_INDEXED_VIEW_METHOD_CONST
150 {
151  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
152  return IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,IvcIndex,typename IvcType<Indices>::type>
153  (derived(), IvcIndex(0), ivcCol(indices));
154 }
155 
156 template<typename Indices>
157 typename internal::enable_if<
158  (!IsRowMajor) && (!(internal::get_compile_time_incr<typename IvcType<Indices>::type>::value==1 || internal::is_valid_index_type<Indices>::value)),
159  IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,typename IvcType<Indices>::type,IvcIndex> >::type
160 operator()(const Indices& indices) EIGEN_INDEXED_VIEW_METHOD_CONST
161 {
162  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
163  return IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,typename IvcType<Indices>::type,IvcIndex>
164  (derived(), ivcRow(indices), IvcIndex(0));
165 }
166 
167 template<typename Indices>
168 typename internal::enable_if<
169  (internal::get_compile_time_incr<typename IvcType<Indices>::type>::value==1) && (!internal::is_valid_index_type<Indices>::value) && (!Symbolic::is_symbolic<Indices>::value),
170  VectorBlock<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,internal::array_size<Indices>::value> >::type
171 operator()(const Indices& indices) EIGEN_INDEXED_VIEW_METHOD_CONST
172 {
173  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
174  typename IvcType<Indices>::type actualIndices = ivcSize(indices);
175  return VectorBlock<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,internal::array_size<Indices>::value>
176  (derived(), internal::first(actualIndices), internal::size(actualIndices));
177 }
178 
179 template<typename IndexType>
180 typename internal::enable_if<Symbolic::is_symbolic<IndexType>::value, CoeffReturnType >::type
181 operator()(const IndexType& id) EIGEN_INDEXED_VIEW_METHOD_CONST
182 {
183  return Base::operator()(internal::eval_expr_given_size(id,size()));
184 }
185 
186 #if EIGEN_HAS_STATIC_ARRAY_TEMPLATE
187 
188 template<typename IndicesT, std::size_t IndicesN>
189 typename internal::enable_if<IsRowMajor,
190  IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,IvcIndex,const IndicesT (&)[IndicesN]> >::type
191 operator()(const IndicesT (&indices)[IndicesN]) EIGEN_INDEXED_VIEW_METHOD_CONST
192 {
193  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
194  return IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,IvcIndex,const IndicesT (&)[IndicesN]>
195  (derived(), IvcIndex(0), indices);
196 }
197 
198 template<typename IndicesT, std::size_t IndicesN>
199 typename internal::enable_if<!IsRowMajor,
200  IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,const IndicesT (&)[IndicesN],IvcIndex> >::type
201 operator()(const IndicesT (&indices)[IndicesN]) EIGEN_INDEXED_VIEW_METHOD_CONST
202 {
203  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
204  return IndexedView<EIGEN_INDEXED_VIEW_METHOD_CONST Derived,const IndicesT (&)[IndicesN],IvcIndex>
205  (derived(), indices, IvcIndex(0));
206 }
207 
208 #endif // EIGEN_HAS_STATIC_ARRAY_TEMPLATE
209 
210 #undef EIGEN_INDEXED_VIEW_METHOD_CONST
211 #undef EIGEN_INDEXED_VIEW_METHOD_TYPE
212 
213 #ifndef EIGEN_INDEXED_VIEW_METHOD_2ND_PASS
214 #define EIGEN_INDEXED_VIEW_METHOD_2ND_PASS
215 #include "IndexedViewMethods.h"
216 #undef EIGEN_INDEXED_VIEW_METHOD_2ND_PASS
217 #endif
218 
219 #else // EIGEN_PARSED_BY_DOXYGEN
220 
255 template<typename RowIndices, typename ColIndices>
256 IndexedView_or_Block
257 operator()(const RowIndices& rowIndices, const ColIndices& colIndices);
258 
263 template<typename Indices>
264 IndexedView_or_VectorBlock
265 operator()(const Indices& indices);
266 
267 #endif // EIGEN_PARSED_BY_DOXYGEN