Eigen  3.4.90 (git rev 67eeba6e720c5745abc77ae6c92ce0a44aa7b7ae)
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 #include "./InternalHeaderCheck.h"
14 
15 namespace Eigen {
16 
17 namespace internal {
18 
19 // helper function for dot(). The problem is that if we put that in the body of dot(), then upon calling dot
20 // with mismatched types, the compiler emits errors about failing to instantiate cwiseProduct BEFORE
21 // looking at the static assertions. Thus this is a trick to get better compile errors.
22 template<typename T, typename U,
23  bool NeedToTranspose = T::IsVectorAtCompileTime && U::IsVectorAtCompileTime &&
24  ((int(T::RowsAtCompileTime) == 1 && int(U::ColsAtCompileTime) == 1) ||
25  (int(T::ColsAtCompileTime) == 1 && int(U::RowsAtCompileTime) == 1))>
26 struct dot_nocheck
27 {
28  typedef scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> conj_prod;
29  typedef typename conj_prod::result_type ResScalar;
30  EIGEN_DEVICE_FUNC
31  EIGEN_STRONG_INLINE
32  static ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b)
33  {
34  return a.template binaryExpr<conj_prod>(b).sum();
35  }
36 };
37 
38 template<typename T, typename U>
39 struct dot_nocheck<T, U, true>
40 {
41  typedef scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> conj_prod;
42  typedef typename conj_prod::result_type ResScalar;
43  EIGEN_DEVICE_FUNC
44  EIGEN_STRONG_INLINE
45  static ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b)
46  {
47  return a.transpose().template binaryExpr<conj_prod>(b).sum();
48  }
49 };
50 
51 } // end namespace internal
52 
64 template<typename Derived>
65 template<typename OtherDerived>
66 EIGEN_DEVICE_FUNC
67 EIGEN_STRONG_INLINE
68 typename ScalarBinaryOpTraits<typename internal::traits<Derived>::Scalar,typename internal::traits<OtherDerived>::Scalar>::ReturnType
70 {
71  EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived)
72  EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived)
73  EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
74 #if !(defined(EIGEN_NO_STATIC_ASSERT) && defined(EIGEN_NO_DEBUG))
75  typedef internal::scalar_conj_product_op<Scalar,typename OtherDerived::Scalar> func;
76  EIGEN_CHECK_BINARY_COMPATIBILIY(func,Scalar,typename OtherDerived::Scalar);
77 #endif
78 
79  eigen_assert(size() == other.size());
80 
81  return internal::dot_nocheck<Derived,OtherDerived>::run(*this, other);
82 }
83 
84 //---------- implementation of L2 norm and related functions ----------
85 
92 template<typename Derived>
93 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename NumTraits<typename internal::traits<Derived>::Scalar>::Real MatrixBase<Derived>::squaredNorm() const
94 {
95  return numext::real((*this).cwiseAbs2().sum());
96 }
97 
104 template<typename Derived>
105 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename NumTraits<typename internal::traits<Derived>::Scalar>::Real MatrixBase<Derived>::norm() const
106 {
107  return numext::sqrt(squaredNorm());
108 }
109 
119 template<typename Derived>
120 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const typename MatrixBase<Derived>::PlainObject
122 {
123  typedef typename internal::nested_eval<Derived,2>::type Nested_;
124  Nested_ n(derived());
125  RealScalar z = n.squaredNorm();
126  // NOTE: after extensive benchmarking, this conditional does not impact performance, at least on recent x86 CPU
127  if(z>RealScalar(0))
128  return n / numext::sqrt(z);
129  else
130  return n;
131 }
132 
141 template<typename Derived>
142 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void MatrixBase<Derived>::normalize()
143 {
144  RealScalar z = squaredNorm();
145  // NOTE: after extensive benchmarking, this conditional does not impact performance, at least on recent x86 CPU
146  if(z>RealScalar(0))
147  derived() /= numext::sqrt(z);
148 }
149 
162 template<typename Derived>
163 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const typename MatrixBase<Derived>::PlainObject
165 {
166  typedef typename internal::nested_eval<Derived,3>::type Nested_;
167  Nested_ n(derived());
168  RealScalar w = n.cwiseAbs().maxCoeff();
169  RealScalar z = (n/w).squaredNorm();
170  if(z>RealScalar(0))
171  return n / (numext::sqrt(z)*w);
172  else
173  return n;
174 }
175 
187 template<typename Derived>
188 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void MatrixBase<Derived>::stableNormalize()
189 {
190  RealScalar w = cwiseAbs().maxCoeff();
191  RealScalar z = (derived()/w).squaredNorm();
192  if(z>RealScalar(0))
193  derived() /= numext::sqrt(z)*w;
194 }
195 
196 //---------- implementation of other norms ----------
197 
198 namespace internal {
199 
200 template<typename Derived, int p>
201 struct lpNorm_selector
202 {
203  typedef typename NumTraits<typename traits<Derived>::Scalar>::Real RealScalar;
204  EIGEN_DEVICE_FUNC
205  static inline RealScalar run(const MatrixBase<Derived>& m)
206  {
207  EIGEN_USING_STD(pow)
208  return pow(m.cwiseAbs().array().pow(p).sum(), RealScalar(1)/p);
209  }
210 };
211 
212 template<typename Derived>
213 struct lpNorm_selector<Derived, 1>
214 {
215  EIGEN_DEVICE_FUNC
216  static inline typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m)
217  {
218  return m.cwiseAbs().sum();
219  }
220 };
221 
222 template<typename Derived>
223 struct lpNorm_selector<Derived, 2>
224 {
225  EIGEN_DEVICE_FUNC
226  static inline typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m)
227  {
228  return m.norm();
229  }
230 };
231 
232 template<typename Derived>
233 struct lpNorm_selector<Derived, Infinity>
234 {
235  typedef typename NumTraits<typename traits<Derived>::Scalar>::Real RealScalar;
236  EIGEN_DEVICE_FUNC
237  static inline RealScalar run(const MatrixBase<Derived>& m)
238  {
239  if(Derived::SizeAtCompileTime==0 || (Derived::SizeAtCompileTime==Dynamic && m.size()==0))
240  return RealScalar(0);
241  return m.cwiseAbs().maxCoeff();
242  }
243 };
244 
245 } // end namespace internal
246 
257 template<typename Derived>
258 template<int p>
259 #ifndef EIGEN_PARSED_BY_DOXYGEN
260 EIGEN_DEVICE_FUNC inline typename NumTraits<typename internal::traits<Derived>::Scalar>::Real
261 #else
262 EIGEN_DEVICE_FUNC MatrixBase<Derived>::RealScalar
263 #endif
265 {
266  return internal::lpNorm_selector<Derived, p>::run(*this);
267 }
268 
269 //---------- implementation of isOrthogonal / isUnitary ----------
270 
277 template<typename Derived>
278 template<typename OtherDerived>
280 (const MatrixBase<OtherDerived>& other, const RealScalar& prec) const
281 {
282  typename internal::nested_eval<Derived,2>::type nested(derived());
283  typename internal::nested_eval<OtherDerived,2>::type otherNested(other.derived());
284  return numext::abs2(nested.dot(otherNested)) <= prec * prec * nested.squaredNorm() * otherNested.squaredNorm();
285 }
286 
298 template<typename Derived>
299 bool MatrixBase<Derived>::isUnitary(const RealScalar& prec) const
300 {
301  typename internal::nested_eval<Derived,1>::type self(derived());
302  for(Index i = 0; i < cols(); ++i)
303  {
304  if(!internal::isApprox(self.col(i).squaredNorm(), static_cast<RealScalar>(1), prec))
305  return false;
306  for(Index j = 0; j < i; ++j)
307  if(!internal::isMuchSmallerThan(self.col(i).dot(self.col(j)), static_cast<Scalar>(1), prec))
308  return false;
309  }
310  return true;
311 }
312 
313 } // end namespace Eigen
314 
315 #endif // EIGEN_DOT_H
internal::traits< Derived >::Scalar Scalar
Definition: DenseBase.h:61
Scalar sum() const
Definition: Redux.h:461
Derived & derived()
Definition: EigenBase.h:48
EIGEN_CONSTEXPR Index size() const EIGEN_NOEXCEPT
Definition: EigenBase.h:69
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:52
void stableNormalize()
Definition: Dot.h:188
RealScalar norm() const
Definition: Dot.h:105
ArrayWrapper< Derived > array()
Definition: MatrixBase.h:315
const PlainObject stableNormalized() const
Definition: Dot.h:164
const PlainObject normalized() const
Definition: Dot.h:121
RealScalar lpNorm() const
Definition: Dot.h:264
bool isUnitary(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Dot.h:299
RealScalar squaredNorm() const
Definition: Dot.h:93
void normalize()
Definition: Dot.h:142
ScalarBinaryOpTraits< typename internal::traits< Derived >::Scalar, typename internal::traits< OtherDerived >::Scalar >::ReturnType dot(const MatrixBase< OtherDerived > &other) const
Definition: Dot.h:69
bool isOrthogonal(const MatrixBase< OtherDerived > &other, const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Dot.h:280
Namespace containing all symbols from the Eigen library.
Definition: Core:139
const int Infinity
Definition: Constants.h:38
const int Dynamic
Definition: Constants.h:24
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:41
Holds information about the various numeric (i.e. scalar) types allowed by Eigen.
Definition: NumTraits.h:231