Eigen  3.2.93
Dot.h
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 // 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_DOT_H
11 #define EIGEN_DOT_H
12 
13 namespace Eigen {
14 
15 namespace internal {
16 
17 // helper function for dot(). The problem is that if we put that in the body of dot(), then upon calling dot
18 // with mismatched types, the compiler emits errors about failing to instantiate cwiseProduct BEFORE
19 // looking at the static assertions. Thus this is a trick to get better compile errors.
20 template<typename T, typename U,
21 // the NeedToTranspose condition here is taken straight from Assign.h
22  bool NeedToTranspose = T::IsVectorAtCompileTime
23  && U::IsVectorAtCompileTime
24  && ((int(T::RowsAtCompileTime) == 1 && int(U::ColsAtCompileTime) == 1)
25  | // FIXME | instead of || to please GCC 4.4.0 stupid warning "suggest parentheses around &&".
26  // revert to || as soon as not needed anymore.
27  (int(T::ColsAtCompileTime) == 1 && int(U::RowsAtCompileTime) == 1))
28 >
29 struct dot_nocheck
30 {
31  typedef scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> conj_prod;
32  typedef typename conj_prod::result_type ResScalar;
33  EIGEN_DEVICE_FUNC
34  static inline ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b)
35  {
36  return a.template binaryExpr<conj_prod>(b).sum();
37  }
38 };
39 
40 template<typename T, typename U>
41 struct dot_nocheck<T, U, true>
42 {
43  typedef scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> conj_prod;
44  typedef typename conj_prod::result_type ResScalar;
45  EIGEN_DEVICE_FUNC
46  static inline ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b)
47  {
48  return a.transpose().template binaryExpr<conj_prod>(b).sum();
49  }
50 };
51 
52 } // end namespace internal
53 
64 template<typename Derived>
65 template<typename OtherDerived>
66 EIGEN_DEVICE_FUNC
67 typename ScalarBinaryOpTraits<typename internal::traits<Derived>::Scalar,typename internal::traits<OtherDerived>::Scalar>::ReturnType
69 {
70  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
71  EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived)
72  EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
73  typedef internal::scalar_conj_product_op<Scalar,typename OtherDerived::Scalar> func;
74  EIGEN_CHECK_BINARY_COMPATIBILIY(func,Scalar,typename OtherDerived::Scalar);
75 
76  eigen_assert(size() == other.size());
77 
78  return internal::dot_nocheck<Derived,OtherDerived>::run(*this, other);
79 }
80 
81 //---------- implementation of L2 norm and related functions ----------
82 
89 template<typename Derived>
91 {
92  return numext::real((*this).cwiseAbs2().sum());
93 }
94 
101 template<typename Derived>
103 {
104  return numext::sqrt(squaredNorm());
105 }
106 
116 template<typename Derived>
117 inline const typename MatrixBase<Derived>::PlainObject
119 {
120  typedef typename internal::nested_eval<Derived,2>::type _Nested;
121  _Nested n(derived());
122  RealScalar z = n.squaredNorm();
123  // NOTE: after extensive benchmarking, this conditional does not impact performance, at least on recent x86 CPU
124  if(z>RealScalar(0))
125  return n / numext::sqrt(z);
126  else
127  return n;
128 }
129 
138 template<typename Derived>
140 {
141  RealScalar z = squaredNorm();
142  // NOTE: after extensive benchmarking, this conditional does not impact performance, at least on recent x86 CPU
143  if(z>RealScalar(0))
144  derived() /= numext::sqrt(z);
145 }
146 
159 template<typename Derived>
160 inline const typename MatrixBase<Derived>::PlainObject
162 {
163  typedef typename internal::nested_eval<Derived,3>::type _Nested;
164  _Nested n(derived());
165  RealScalar w = n.cwiseAbs().maxCoeff();
166  RealScalar z = (n/w).squaredNorm();
167  if(z>RealScalar(0))
168  return n / (numext::sqrt(z)*w);
169  else
170  return n;
171 }
172 
184 template<typename Derived>
186 {
187  RealScalar w = cwiseAbs().maxCoeff();
188  RealScalar z = (derived()/w).squaredNorm();
189  if(z>RealScalar(0))
190  derived() /= numext::sqrt(z)*w;
191 }
192 
193 //---------- implementation of other norms ----------
194 
195 namespace internal {
196 
197 template<typename Derived, int p>
198 struct lpNorm_selector
199 {
200  typedef typename NumTraits<typename traits<Derived>::Scalar>::Real RealScalar;
201  EIGEN_DEVICE_FUNC
202  static inline RealScalar run(const MatrixBase<Derived>& m)
203  {
204  EIGEN_USING_STD_MATH(pow)
205  return pow(m.cwiseAbs().array().pow(p).sum(), RealScalar(1)/p);
206  }
207 };
208 
209 template<typename Derived>
210 struct lpNorm_selector<Derived, 1>
211 {
212  EIGEN_DEVICE_FUNC
213  static inline typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m)
214  {
215  return m.cwiseAbs().sum();
216  }
217 };
218 
219 template<typename Derived>
220 struct lpNorm_selector<Derived, 2>
221 {
222  EIGEN_DEVICE_FUNC
223  static inline typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m)
224  {
225  return m.norm();
226  }
227 };
228 
229 template<typename Derived>
230 struct lpNorm_selector<Derived, Infinity>
231 {
232  typedef typename NumTraits<typename traits<Derived>::Scalar>::Real RealScalar;
233  EIGEN_DEVICE_FUNC
234  static inline RealScalar run(const MatrixBase<Derived>& m)
235  {
236  if(Derived::SizeAtCompileTime==0 || (Derived::SizeAtCompileTime==Dynamic && m.size()==0))
237  return RealScalar(0);
238  return m.cwiseAbs().maxCoeff();
239  }
240 };
241 
242 } // end namespace internal
243 
254 template<typename Derived>
255 template<int p>
256 #ifndef EIGEN_PARSED_BY_DOXYGEN
258 #else
259 MatrixBase<Derived>::RealScalar
260 #endif
262 {
263  return internal::lpNorm_selector<Derived, p>::run(*this);
264 }
265 
266 //---------- implementation of isOrthogonal / isUnitary ----------
267 
274 template<typename Derived>
275 template<typename OtherDerived>
277 (const MatrixBase<OtherDerived>& other, const RealScalar& prec) const
278 {
279  typename internal::nested_eval<Derived,2>::type nested(derived());
280  typename internal::nested_eval<OtherDerived,2>::type otherNested(other.derived());
281  return numext::abs2(nested.dot(otherNested)) <= prec * prec * nested.squaredNorm() * otherNested.squaredNorm();
282 }
283 
295 template<typename Derived>
296 bool MatrixBase<Derived>::isUnitary(const RealScalar& prec) const
297 {
298  typename internal::nested_eval<Derived,1>::type self(derived());
299  for(Index i = 0; i < cols(); ++i)
300  {
301  if(!internal::isApprox(self.col(i).squaredNorm(), static_cast<RealScalar>(1), prec))
302  return false;
303  for(Index j = 0; j < i; ++j)
304  if(!internal::isMuchSmallerThan(self.col(i).dot(self.col(j)), static_cast<Scalar>(1), prec))
305  return false;
306  }
307  return true;
308 }
309 
310 } // end namespace Eigen
311 
312 #endif // EIGEN_DOT_H
internal::traits< Derived >::Scalar Scalar
Definition: DenseBase.h:66
bool isUnitary(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Dot.h:296
Namespace containing all symbols from the Eigen library.
Definition: Core:271
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:167
const PlainObject normalized() const
Definition: Dot.h:118
Derived & derived()
Definition: EigenBase.h:44
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:37
void normalize()
Definition: Dot.h:139
RealScalar norm() const
Definition: Dot.h:102
void stableNormalize()
Definition: Dot.h:185
RealScalar squaredNorm() const
Definition: Dot.h:90
Index size() const
Definition: EigenBase.h:65
const PlainObject stableNormalized() const
Definition: Dot.h:161
Definition: Eigen_Colamd.h:50
ScalarBinaryOpTraits< typename internal::traits< Derived >::Scalar, typename internal::traits< OtherDerived >::Scalar >::ReturnType dot(const MatrixBase< OtherDerived > &other) const
Definition: Dot.h:68
RealScalar lpNorm() const
Definition: Dot.h:261
bool isOrthogonal(const MatrixBase< OtherDerived > &other, const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Dot.h:277
const int Dynamic
Definition: Constants.h:21
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const CwiseAbsReturnType cwiseAbs() const
Definition: MatrixBase.h:30