Another set of UI messages fixes and tweaks! No functional changes.
[blender.git] / extern / Eigen3 / Eigen / src / Core / Visitor.h
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 // Eigen is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 3 of the License, or (at your option) any later version.
10 //
11 // Alternatively, you can redistribute it and/or
12 // modify it under the terms of the GNU General Public License as
13 // published by the Free Software Foundation; either version 2 of
14 // the License, or (at your option) any later version.
15 //
16 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
17 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU Lesser General Public
22 // License and a copy of the GNU General Public License along with
23 // Eigen. If not, see <http://www.gnu.org/licenses/>.
24
25 #ifndef EIGEN_VISITOR_H
26 #define EIGEN_VISITOR_H
27
28 namespace internal {
29
30 template<typename Visitor, typename Derived, int UnrollCount>
31 struct visitor_impl
32 {
33   enum {
34     col = (UnrollCount-1) / Derived::RowsAtCompileTime,
35     row = (UnrollCount-1) % Derived::RowsAtCompileTime
36   };
37
38   inline static void run(const Derived &mat, Visitor& visitor)
39   {
40     visitor_impl<Visitor, Derived, UnrollCount-1>::run(mat, visitor);
41     visitor(mat.coeff(row, col), row, col);
42   }
43 };
44
45 template<typename Visitor, typename Derived>
46 struct visitor_impl<Visitor, Derived, 1>
47 {
48   inline static void run(const Derived &mat, Visitor& visitor)
49   {
50     return visitor.init(mat.coeff(0, 0), 0, 0);
51   }
52 };
53
54 template<typename Visitor, typename Derived>
55 struct visitor_impl<Visitor, Derived, Dynamic>
56 {
57   typedef typename Derived::Index Index;
58   inline static void run(const Derived& mat, Visitor& visitor)
59   {
60     visitor.init(mat.coeff(0,0), 0, 0);
61     for(Index i = 1; i < mat.rows(); ++i)
62       visitor(mat.coeff(i, 0), i, 0);
63     for(Index j = 1; j < mat.cols(); ++j)
64       for(Index i = 0; i < mat.rows(); ++i)
65         visitor(mat.coeff(i, j), i, j);
66   }
67 };
68
69 } // end namespace internal
70
71 /** Applies the visitor \a visitor to the whole coefficients of the matrix or vector.
72   *
73   * The template parameter \a Visitor is the type of the visitor and provides the following interface:
74   * \code
75   * struct MyVisitor {
76   *   // called for the first coefficient
77   *   void init(const Scalar& value, Index i, Index j);
78   *   // called for all other coefficients
79   *   void operator() (const Scalar& value, Index i, Index j);
80   * };
81   * \endcode
82   *
83   * \note compared to one or two \em for \em loops, visitors offer automatic
84   * unrolling for small fixed size matrix.
85   *
86   * \sa minCoeff(Index*,Index*), maxCoeff(Index*,Index*), DenseBase::redux()
87   */
88 template<typename Derived>
89 template<typename Visitor>
90 void DenseBase<Derived>::visit(Visitor& visitor) const
91 {
92   enum { unroll = SizeAtCompileTime != Dynamic
93                    && CoeffReadCost != Dynamic
94                    && (SizeAtCompileTime == 1 || internal::functor_traits<Visitor>::Cost != Dynamic)
95                    && SizeAtCompileTime * CoeffReadCost + (SizeAtCompileTime-1) * internal::functor_traits<Visitor>::Cost
96                       <= EIGEN_UNROLLING_LIMIT };
97   return internal::visitor_impl<Visitor, Derived,
98       unroll ? int(SizeAtCompileTime) : Dynamic
99     >::run(derived(), visitor);
100 }
101
102 namespace internal {
103
104 /** \internal
105   * \brief Base class to implement min and max visitors
106   */
107 template <typename Derived>
108 struct coeff_visitor
109 {
110   typedef typename Derived::Index Index;
111   typedef typename Derived::Scalar Scalar;
112   Index row, col;
113   Scalar res;
114   inline void init(const Scalar& value, Index i, Index j)
115   {
116     res = value;
117     row = i;
118     col = j;
119   }
120 };
121
122 /** \internal
123   * \brief Visitor computing the min coefficient with its value and coordinates
124   *
125   * \sa DenseBase::minCoeff(Index*, Index*)
126   */
127 template <typename Derived>
128 struct min_coeff_visitor : coeff_visitor<Derived>
129 {
130   typedef typename Derived::Index Index;
131   typedef typename Derived::Scalar Scalar;
132   void operator() (const Scalar& value, Index i, Index j)
133   {
134     if(value < this->res)
135     {
136       this->res = value;
137       this->row = i;
138       this->col = j;
139     }
140   }
141 };
142
143 template<typename Scalar>
144 struct functor_traits<min_coeff_visitor<Scalar> > {
145   enum {
146     Cost = NumTraits<Scalar>::AddCost
147   };
148 };
149
150 /** \internal
151   * \brief Visitor computing the max coefficient with its value and coordinates
152   *
153   * \sa DenseBase::maxCoeff(Index*, Index*)
154   */
155 template <typename Derived>
156 struct max_coeff_visitor : coeff_visitor<Derived>
157 {
158   typedef typename Derived::Index Index;
159   typedef typename Derived::Scalar Scalar;
160   void operator() (const Scalar& value, Index i, Index j)
161   {
162     if(value > this->res)
163     {
164       this->res = value;
165       this->row = i;
166       this->col = j;
167     }
168   }
169 };
170
171 template<typename Scalar>
172 struct functor_traits<max_coeff_visitor<Scalar> > {
173   enum {
174     Cost = NumTraits<Scalar>::AddCost
175   };
176 };
177
178 } // end namespace internal
179
180 /** \returns the minimum of all coefficients of *this
181   * and puts in *row and *col its location.
182   *
183   * \sa DenseBase::minCoeff(Index*), DenseBase::maxCoeff(Index*,Index*), DenseBase::visitor(), DenseBase::minCoeff()
184   */
185 template<typename Derived>
186 template<typename IndexType>
187 typename internal::traits<Derived>::Scalar
188 DenseBase<Derived>::minCoeff(IndexType* row, IndexType* col) const
189 {
190   internal::min_coeff_visitor<Derived> minVisitor;
191   this->visit(minVisitor);
192   *row = minVisitor.row;
193   if (col) *col = minVisitor.col;
194   return minVisitor.res;
195 }
196
197 /** \returns the minimum of all coefficients of *this
198   * and puts in *index its location.
199   *
200   * \sa DenseBase::minCoeff(IndexType*,IndexType*), DenseBase::maxCoeff(IndexType*,IndexType*), DenseBase::visitor(), DenseBase::minCoeff()
201   */
202 template<typename Derived>
203 template<typename IndexType>
204 typename internal::traits<Derived>::Scalar
205 DenseBase<Derived>::minCoeff(IndexType* index) const
206 {
207   EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
208   internal::min_coeff_visitor<Derived> minVisitor;
209   this->visit(minVisitor);
210   *index = (RowsAtCompileTime==1) ? minVisitor.col : minVisitor.row;
211   return minVisitor.res;
212 }
213
214 /** \returns the maximum of all coefficients of *this
215   * and puts in *row and *col its location.
216   *
217   * \sa DenseBase::minCoeff(IndexType*,IndexType*), DenseBase::visitor(), DenseBase::maxCoeff()
218   */
219 template<typename Derived>
220 template<typename IndexType>
221 typename internal::traits<Derived>::Scalar
222 DenseBase<Derived>::maxCoeff(IndexType* row, IndexType* col) const
223 {
224   internal::max_coeff_visitor<Derived> maxVisitor;
225   this->visit(maxVisitor);
226   *row = maxVisitor.row;
227   if (col) *col = maxVisitor.col;
228   return maxVisitor.res;
229 }
230
231 /** \returns the maximum of all coefficients of *this
232   * and puts in *index its location.
233   *
234   * \sa DenseBase::maxCoeff(IndexType*,IndexType*), DenseBase::minCoeff(IndexType*,IndexType*), DenseBase::visitor(), DenseBase::maxCoeff()
235   */
236 template<typename Derived>
237 template<typename IndexType>
238 typename internal::traits<Derived>::Scalar
239 DenseBase<Derived>::maxCoeff(IndexType* index) const
240 {
241   EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
242   internal::max_coeff_visitor<Derived> maxVisitor;
243   this->visit(maxVisitor);
244   *index = (RowsAtCompileTime==1) ? maxVisitor.col : maxVisitor.row;
245   return maxVisitor.res;
246 }
247
248 #endif // EIGEN_VISITOR_H