~siretart/ubuntu/utopic/blender/libav10

« back to all changes in this revision

Viewing changes to extern/Eigen3/Eigen/src/Core/Dot.h

  • Committer: Package Import Robot
  • Author(s): Matteo F. Vescovi
  • Date: 2012-07-23 08:54:18 UTC
  • mfrom: (14.2.16 sid)
  • mto: (14.2.19 sid)
  • mto: This revision was merged to the branch mainline in revision 42.
  • Revision ID: package-import@ubuntu.com-20120723085418-9foz30v6afaf5ffs
Tags: 2.63a-2
* debian/: Cycles support added (Closes: #658075)
  For now, this top feature has been enabled only
  on [any-amd64 any-i386] architectures because
  of OpenImageIO failing on all others
* debian/: scripts installation path changed
  from /usr/lib to /usr/share:
  + debian/patches/: patchset re-worked for path changing
  + debian/control: "Breaks" field added on yafaray-exporter

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
// This file is part of Eigen, a lightweight C++ template library
 
2
// for linear algebra.
 
3
//
 
4
// Copyright (C) 2006-2008, 2010 Benoit Jacob <jacob.benoit.1@gmail.com>
 
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_DOT_H
 
26
#define EIGEN_DOT_H
 
27
 
 
28
namespace internal {
 
29
 
 
30
// helper function for dot(). The problem is that if we put that in the body of dot(), then upon calling dot
 
31
// with mismatched types, the compiler emits errors about failing to instantiate cwiseProduct BEFORE
 
32
// looking at the static assertions. Thus this is a trick to get better compile errors.
 
33
template<typename T, typename U,
 
34
// the NeedToTranspose condition here is taken straight from Assign.h
 
35
         bool NeedToTranspose = T::IsVectorAtCompileTime
 
36
                && U::IsVectorAtCompileTime
 
37
                && ((int(T::RowsAtCompileTime) == 1 && int(U::ColsAtCompileTime) == 1)
 
38
                      |  // FIXME | instead of || to please GCC 4.4.0 stupid warning "suggest parentheses around &&".
 
39
                         // revert to || as soon as not needed anymore.
 
40
                    (int(T::ColsAtCompileTime) == 1 && int(U::RowsAtCompileTime) == 1))
 
41
>
 
42
struct dot_nocheck
 
43
{
 
44
  typedef typename scalar_product_traits<typename traits<T>::Scalar,typename traits<U>::Scalar>::ReturnType ResScalar;
 
45
  static inline ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b)
 
46
  {
 
47
    return a.template binaryExpr<scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> >(b).sum();
 
48
  }
 
49
};
 
50
 
 
51
template<typename T, typename U>
 
52
struct dot_nocheck<T, U, true>
 
53
{
 
54
  typedef typename scalar_product_traits<typename traits<T>::Scalar,typename traits<U>::Scalar>::ReturnType ResScalar;
 
55
  static inline ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b)
 
56
  {
 
57
    return a.transpose().template binaryExpr<scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> >(b).sum();
 
58
  }
 
59
};
 
60
 
 
61
} // end namespace internal
 
62
 
 
63
/** \returns the dot product of *this with other.
 
64
  *
 
65
  * \only_for_vectors
 
66
  *
 
67
  * \note If the scalar type is complex numbers, then this function returns the hermitian
 
68
  * (sesquilinear) dot product, conjugate-linear in the first variable and linear in the
 
69
  * second variable.
 
70
  *
 
71
  * \sa squaredNorm(), norm()
 
72
  */
 
73
template<typename Derived>
 
74
template<typename OtherDerived>
 
75
typename internal::scalar_product_traits<typename internal::traits<Derived>::Scalar,typename internal::traits<OtherDerived>::Scalar>::ReturnType
 
76
MatrixBase<Derived>::dot(const MatrixBase<OtherDerived>& other) const
 
77
{
 
78
  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
 
79
  EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived)
 
80
  EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
 
81
  typedef internal::scalar_conj_product_op<Scalar,typename OtherDerived::Scalar> func;
 
82
  EIGEN_CHECK_BINARY_COMPATIBILIY(func,Scalar,typename OtherDerived::Scalar);
 
83
 
 
84
  eigen_assert(size() == other.size());
 
85
 
 
86
  return internal::dot_nocheck<Derived,OtherDerived>::run(*this, other);
 
87
}
 
88
 
 
89
#ifdef EIGEN2_SUPPORT
 
90
/** \returns the dot product of *this with other, with the Eigen2 convention that the dot product is linear in the first variable
 
91
  * (conjugating the second variable). Of course this only makes a difference in the complex case.
 
92
  *
 
93
  * This method is only available in EIGEN2_SUPPORT mode.
 
94
  *
 
95
  * \only_for_vectors
 
96
  *
 
97
  * \sa dot()
 
98
  */
 
99
template<typename Derived>
 
100
template<typename OtherDerived>
 
101
typename internal::traits<Derived>::Scalar
 
102
MatrixBase<Derived>::eigen2_dot(const MatrixBase<OtherDerived>& other) const
 
103
{
 
104
  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
 
105
  EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived)
 
106
  EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
 
107
  EIGEN_STATIC_ASSERT((internal::is_same<Scalar, typename OtherDerived::Scalar>::value),
 
108
    YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
 
109
 
 
110
  eigen_assert(size() == other.size());
 
111
 
 
112
  return internal::dot_nocheck<OtherDerived,Derived>::run(other,*this);
 
113
}
 
114
#endif
 
115
 
 
116
 
 
117
//---------- implementation of L2 norm and related functions ----------
 
118
 
 
119
/** \returns, for vectors, the squared \em l2 norm of \c *this, and for matrices the Frobenius norm.
 
120
  * In both cases, it consists in the sum of the square of all the matrix entries.
 
121
  * For vectors, this is also equals to the dot product of \c *this with itself.
 
122
  *
 
123
  * \sa dot(), norm()
 
124
  */
 
125
template<typename Derived>
 
126
EIGEN_STRONG_INLINE typename NumTraits<typename internal::traits<Derived>::Scalar>::Real MatrixBase<Derived>::squaredNorm() const
 
127
{
 
128
  return internal::real((*this).cwiseAbs2().sum());
 
129
}
 
130
 
 
131
/** \returns, for vectors, the \em l2 norm of \c *this, and for matrices the Frobenius norm.
 
132
  * In both cases, it consists in the square root of the sum of the square of all the matrix entries.
 
133
  * For vectors, this is also equals to the square root of the dot product of \c *this with itself.
 
134
  *
 
135
  * \sa dot(), squaredNorm()
 
136
  */
 
137
template<typename Derived>
 
138
inline typename NumTraits<typename internal::traits<Derived>::Scalar>::Real MatrixBase<Derived>::norm() const
 
139
{
 
140
  return internal::sqrt(squaredNorm());
 
141
}
 
142
 
 
143
/** \returns an expression of the quotient of *this by its own norm.
 
144
  *
 
145
  * \only_for_vectors
 
146
  *
 
147
  * \sa norm(), normalize()
 
148
  */
 
149
template<typename Derived>
 
150
inline const typename MatrixBase<Derived>::PlainObject
 
151
MatrixBase<Derived>::normalized() const
 
152
{
 
153
  typedef typename internal::nested<Derived>::type Nested;
 
154
  typedef typename internal::remove_reference<Nested>::type _Nested;
 
155
  _Nested n(derived());
 
156
  return n / n.norm();
 
157
}
 
158
 
 
159
/** Normalizes the vector, i.e. divides it by its own norm.
 
160
  *
 
161
  * \only_for_vectors
 
162
  *
 
163
  * \sa norm(), normalized()
 
164
  */
 
165
template<typename Derived>
 
166
inline void MatrixBase<Derived>::normalize()
 
167
{
 
168
  *this /= norm();
 
169
}
 
170
 
 
171
//---------- implementation of other norms ----------
 
172
 
 
173
namespace internal {
 
174
 
 
175
template<typename Derived, int p>
 
176
struct lpNorm_selector
 
177
{
 
178
  typedef typename NumTraits<typename traits<Derived>::Scalar>::Real RealScalar;
 
179
  inline static RealScalar run(const MatrixBase<Derived>& m)
 
180
  {
 
181
    return pow(m.cwiseAbs().array().pow(p).sum(), RealScalar(1)/p);
 
182
  }
 
183
};
 
184
 
 
185
template<typename Derived>
 
186
struct lpNorm_selector<Derived, 1>
 
187
{
 
188
  inline static typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m)
 
189
  {
 
190
    return m.cwiseAbs().sum();
 
191
  }
 
192
};
 
193
 
 
194
template<typename Derived>
 
195
struct lpNorm_selector<Derived, 2>
 
196
{
 
197
  inline static typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m)
 
198
  {
 
199
    return m.norm();
 
200
  }
 
201
};
 
202
 
 
203
template<typename Derived>
 
204
struct lpNorm_selector<Derived, Infinity>
 
205
{
 
206
  inline static typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m)
 
207
  {
 
208
    return m.cwiseAbs().maxCoeff();
 
209
  }
 
210
};
 
211
 
 
212
} // end namespace internal
 
213
 
 
214
/** \returns the \f$ \ell^p \f$ norm of *this, that is, returns the p-th root of the sum of the p-th powers of the absolute values
 
215
  *          of the coefficients of *this. If \a p is the special value \a Eigen::Infinity, this function returns the \f$ \ell^\infty \f$
 
216
  *          norm, that is the maximum of the absolute values of the coefficients of *this.
 
217
  *
 
218
  * \sa norm()
 
219
  */
 
220
template<typename Derived>
 
221
template<int p>
 
222
inline typename NumTraits<typename internal::traits<Derived>::Scalar>::Real
 
223
MatrixBase<Derived>::lpNorm() const
 
224
{
 
225
  return internal::lpNorm_selector<Derived, p>::run(*this);
 
226
}
 
227
 
 
228
//---------- implementation of isOrthogonal / isUnitary ----------
 
229
 
 
230
/** \returns true if *this is approximately orthogonal to \a other,
 
231
  *          within the precision given by \a prec.
 
232
  *
 
233
  * Example: \include MatrixBase_isOrthogonal.cpp
 
234
  * Output: \verbinclude MatrixBase_isOrthogonal.out
 
235
  */
 
236
template<typename Derived>
 
237
template<typename OtherDerived>
 
238
bool MatrixBase<Derived>::isOrthogonal
 
239
(const MatrixBase<OtherDerived>& other, RealScalar prec) const
 
240
{
 
241
  typename internal::nested<Derived,2>::type nested(derived());
 
242
  typename internal::nested<OtherDerived,2>::type otherNested(other.derived());
 
243
  return internal::abs2(nested.dot(otherNested)) <= prec * prec * nested.squaredNorm() * otherNested.squaredNorm();
 
244
}
 
245
 
 
246
/** \returns true if *this is approximately an unitary matrix,
 
247
  *          within the precision given by \a prec. In the case where the \a Scalar
 
248
  *          type is real numbers, a unitary matrix is an orthogonal matrix, whence the name.
 
249
  *
 
250
  * \note This can be used to check whether a family of vectors forms an orthonormal basis.
 
251
  *       Indeed, \c m.isUnitary() returns true if and only if the columns (equivalently, the rows) of m form an
 
252
  *       orthonormal basis.
 
253
  *
 
254
  * Example: \include MatrixBase_isUnitary.cpp
 
255
  * Output: \verbinclude MatrixBase_isUnitary.out
 
256
  */
 
257
template<typename Derived>
 
258
bool MatrixBase<Derived>::isUnitary(RealScalar prec) const
 
259
{
 
260
  typename Derived::Nested nested(derived());
 
261
  for(Index i = 0; i < cols(); ++i)
 
262
  {
 
263
    if(!internal::isApprox(nested.col(i).squaredNorm(), static_cast<RealScalar>(1), prec))
 
264
      return false;
 
265
    for(Index j = 0; j < i; ++j)
 
266
      if(!internal::isMuchSmallerThan(nested.col(i).dot(nested.col(j)), static_cast<Scalar>(1), prec))
 
267
        return false;
 
268
  }
 
269
  return true;
 
270
}
 
271
 
 
272
#endif // EIGEN_DOT_H