• Home
  • History
  • Annotate
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
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 //
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_VISITOR_H
11 #define EIGEN_VISITOR_H
12 
13 namespace Eigen {
14 
15 namespace internal {
16 
17 template<typename Visitor, typename Derived, int UnrollCount>
18 struct visitor_impl
19 {
20   enum {
21     col = (UnrollCount-1) / Derived::RowsAtCompileTime,
22     row = (UnrollCount-1) % Derived::RowsAtCompileTime
23   };
24 
runvisitor_impl25   static inline void run(const Derived &mat, Visitor& visitor)
26   {
27     visitor_impl<Visitor, Derived, UnrollCount-1>::run(mat, visitor);
28     visitor(mat.coeff(row, col), row, col);
29   }
30 };
31 
32 template<typename Visitor, typename Derived>
33 struct visitor_impl<Visitor, Derived, 1>
34 {
35   static inline void run(const Derived &mat, Visitor& visitor)
36   {
37     return visitor.init(mat.coeff(0, 0), 0, 0);
38   }
39 };
40 
41 template<typename Visitor, typename Derived>
42 struct visitor_impl<Visitor, Derived, Dynamic>
43 {
44   typedef typename Derived::Index Index;
45   static inline void run(const Derived& mat, Visitor& visitor)
46   {
47     visitor.init(mat.coeff(0,0), 0, 0);
48     for(Index i = 1; i < mat.rows(); ++i)
49       visitor(mat.coeff(i, 0), i, 0);
50     for(Index j = 1; j < mat.cols(); ++j)
51       for(Index i = 0; i < mat.rows(); ++i)
52         visitor(mat.coeff(i, j), i, j);
53   }
54 };
55 
56 } // end namespace internal
57 
58 /** Applies the visitor \a visitor to the whole coefficients of the matrix or vector.
59   *
60   * The template parameter \a Visitor is the type of the visitor and provides the following interface:
61   * \code
62   * struct MyVisitor {
63   *   // called for the first coefficient
64   *   void init(const Scalar& value, Index i, Index j);
65   *   // called for all other coefficients
66   *   void operator() (const Scalar& value, Index i, Index j);
67   * };
68   * \endcode
69   *
70   * \note compared to one or two \em for \em loops, visitors offer automatic
71   * unrolling for small fixed size matrix.
72   *
73   * \sa minCoeff(Index*,Index*), maxCoeff(Index*,Index*), DenseBase::redux()
74   */
75 template<typename Derived>
76 template<typename Visitor>
77 void DenseBase<Derived>::visit(Visitor& visitor) const
78 {
79   enum { unroll = SizeAtCompileTime != Dynamic
80                    && CoeffReadCost != Dynamic
81                    && (SizeAtCompileTime == 1 || internal::functor_traits<Visitor>::Cost != Dynamic)
82                    && SizeAtCompileTime * CoeffReadCost + (SizeAtCompileTime-1) * internal::functor_traits<Visitor>::Cost
83                       <= EIGEN_UNROLLING_LIMIT };
84   return internal::visitor_impl<Visitor, Derived,
85       unroll ? int(SizeAtCompileTime) : Dynamic
86     >::run(derived(), visitor);
87 }
88 
89 namespace internal {
90 
91 /** \internal
92   * \brief Base class to implement min and max visitors
93   */
94 template <typename Derived>
95 struct coeff_visitor
96 {
97   typedef typename Derived::Index Index;
98   typedef typename Derived::Scalar Scalar;
99   Index row, col;
100   Scalar res;
101   inline void init(const Scalar& value, Index i, Index j)
102   {
103     res = value;
104     row = i;
105     col = j;
106   }
107 };
108 
109 /** \internal
110   * \brief Visitor computing the min coefficient with its value and coordinates
111   *
112   * \sa DenseBase::minCoeff(Index*, Index*)
113   */
114 template <typename Derived>
115 struct min_coeff_visitor : coeff_visitor<Derived>
116 {
117   typedef typename Derived::Index Index;
118   typedef typename Derived::Scalar Scalar;
119   void operator() (const Scalar& value, Index i, Index j)
120   {
121     if(value < this->res)
122     {
123       this->res = value;
124       this->row = i;
125       this->col = j;
126     }
127   }
128 };
129 
130 template<typename Scalar>
131 struct functor_traits<min_coeff_visitor<Scalar> > {
132   enum {
133     Cost = NumTraits<Scalar>::AddCost
134   };
135 };
136 
137 /** \internal
138   * \brief Visitor computing the max coefficient with its value and coordinates
139   *
140   * \sa DenseBase::maxCoeff(Index*, Index*)
141   */
142 template <typename Derived>
143 struct max_coeff_visitor : coeff_visitor<Derived>
144 {
145   typedef typename Derived::Index Index;
146   typedef typename Derived::Scalar Scalar;
147   void operator() (const Scalar& value, Index i, Index j)
148   {
149     if(value > this->res)
150     {
151       this->res = value;
152       this->row = i;
153       this->col = j;
154     }
155   }
156 };
157 
158 template<typename Scalar>
159 struct functor_traits<max_coeff_visitor<Scalar> > {
160   enum {
161     Cost = NumTraits<Scalar>::AddCost
162   };
163 };
164 
165 } // end namespace internal
166 
167 /** \returns the minimum of all coefficients of *this and puts in *row and *col its location.
168   * \warning the result is undefined if \c *this contains NaN.
169   *
170   * \sa DenseBase::minCoeff(Index*), DenseBase::maxCoeff(Index*,Index*), DenseBase::visitor(), DenseBase::minCoeff()
171   */
172 template<typename Derived>
173 template<typename IndexType>
174 typename internal::traits<Derived>::Scalar
175 DenseBase<Derived>::minCoeff(IndexType* rowId, IndexType* colId) const
176 {
177   internal::min_coeff_visitor<Derived> minVisitor;
178   this->visit(minVisitor);
179   *rowId = minVisitor.row;
180   if (colId) *colId = minVisitor.col;
181   return minVisitor.res;
182 }
183 
184 /** \returns the minimum of all coefficients of *this and puts in *index its location.
185   * \warning the result is undefined if \c *this contains NaN.
186   *
187   * \sa DenseBase::minCoeff(IndexType*,IndexType*), DenseBase::maxCoeff(IndexType*,IndexType*), DenseBase::visitor(), DenseBase::minCoeff()
188   */
189 template<typename Derived>
190 template<typename IndexType>
191 typename internal::traits<Derived>::Scalar
192 DenseBase<Derived>::minCoeff(IndexType* index) const
193 {
194   EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
195   internal::min_coeff_visitor<Derived> minVisitor;
196   this->visit(minVisitor);
197   *index = (RowsAtCompileTime==1) ? minVisitor.col : minVisitor.row;
198   return minVisitor.res;
199 }
200 
201 /** \returns the maximum of all coefficients of *this and puts in *row and *col its location.
202   * \warning the result is undefined if \c *this contains NaN.
203   *
204   * \sa DenseBase::minCoeff(IndexType*,IndexType*), DenseBase::visitor(), DenseBase::maxCoeff()
205   */
206 template<typename Derived>
207 template<typename IndexType>
208 typename internal::traits<Derived>::Scalar
209 DenseBase<Derived>::maxCoeff(IndexType* rowPtr, IndexType* colPtr) const
210 {
211   internal::max_coeff_visitor<Derived> maxVisitor;
212   this->visit(maxVisitor);
213   *rowPtr = maxVisitor.row;
214   if (colPtr) *colPtr = maxVisitor.col;
215   return maxVisitor.res;
216 }
217 
218 /** \returns the maximum of all coefficients of *this and puts in *index its location.
219   * \warning the result is undefined if \c *this contains NaN.
220   *
221   * \sa DenseBase::maxCoeff(IndexType*,IndexType*), DenseBase::minCoeff(IndexType*,IndexType*), DenseBase::visitor(), DenseBase::maxCoeff()
222   */
223 template<typename Derived>
224 template<typename IndexType>
225 typename internal::traits<Derived>::Scalar
226 DenseBase<Derived>::maxCoeff(IndexType* index) const
227 {
228   EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
229   internal::max_coeff_visitor<Derived> maxVisitor;
230   this->visit(maxVisitor);
231   *index = (RowsAtCompileTime==1) ? maxVisitor.col : maxVisitor.row;
232   return maxVisitor.res;
233 }
234 
235 } // end namespace Eigen
236 
237 #endif // EIGEN_VISITOR_H
238