![]() |
Eigen
3.3.3
|
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2006-2008, 2010 Benoit Jacob <jacob.benoit.1@gmail.com> 00005 // 00006 // This Source Code Form is subject to the terms of the Mozilla 00007 // Public License v. 2.0. If a copy of the MPL was not distributed 00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00009 00010 #ifndef EIGEN_DOT_H 00011 #define EIGEN_DOT_H 00012 00013 namespace Eigen { 00014 00015 namespace internal { 00016 00017 // helper function for dot(). The problem is that if we put that in the body of dot(), then upon calling dot 00018 // with mismatched types, the compiler emits errors about failing to instantiate cwiseProduct BEFORE 00019 // looking at the static assertions. Thus this is a trick to get better compile errors. 00020 template<typename T, typename U, 00021 // the NeedToTranspose condition here is taken straight from Assign.h 00022 bool NeedToTranspose = T::IsVectorAtCompileTime 00023 && U::IsVectorAtCompileTime 00024 && ((int(T::RowsAtCompileTime) == 1 && int(U::ColsAtCompileTime) == 1) 00025 | // FIXME | instead of || to please GCC 4.4.0 stupid warning "suggest parentheses around &&". 00026 // revert to || as soon as not needed anymore. 00027 (int(T::ColsAtCompileTime) == 1 && int(U::RowsAtCompileTime) == 1)) 00028 > 00029 struct dot_nocheck 00030 { 00031 typedef scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> conj_prod; 00032 typedef typename conj_prod::result_type ResScalar; 00033 EIGEN_DEVICE_FUNC 00034 static inline ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b) 00035 { 00036 return a.template binaryExpr<conj_prod>(b).sum(); 00037 } 00038 }; 00039 00040 template<typename T, typename U> 00041 struct dot_nocheck<T, U, true> 00042 { 00043 typedef scalar_conj_product_op<typename traits<T>::Scalar,typename traits<U>::Scalar> conj_prod; 00044 typedef typename conj_prod::result_type ResScalar; 00045 EIGEN_DEVICE_FUNC 00046 static inline ResScalar run(const MatrixBase<T>& a, const MatrixBase<U>& b) 00047 { 00048 return a.transpose().template binaryExpr<conj_prod>(b).sum(); 00049 } 00050 }; 00051 00052 } // end namespace internal 00053 00065 template<typename Derived> 00066 template<typename OtherDerived> 00067 EIGEN_DEVICE_FUNC 00068 typename ScalarBinaryOpTraits<typename internal::traits<Derived>::Scalar,typename internal::traits<OtherDerived>::Scalar>::ReturnType 00069 MatrixBase<Derived>::dot(const MatrixBase<OtherDerived>& other) const 00070 { 00071 EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived) 00072 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 00073 EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived) 00074 #if !(defined(EIGEN_NO_STATIC_ASSERT) && defined(EIGEN_NO_DEBUG)) 00075 typedef internal::scalar_conj_product_op<Scalar,typename OtherDerived::Scalar> func; 00076 EIGEN_CHECK_BINARY_COMPATIBILIY(func,Scalar,typename OtherDerived::Scalar); 00077 #endif 00078 00079 eigen_assert(size() == other.size()); 00080 00081 return internal::dot_nocheck<Derived,OtherDerived>::run(*this, other); 00082 } 00083 00084 //---------- implementation of L2 norm and related functions ---------- 00085 00092 template<typename Derived> 00093 EIGEN_STRONG_INLINE typename NumTraits<typename internal::traits<Derived>::Scalar>::Real MatrixBase<Derived>::squaredNorm() const 00094 { 00095 return numext::real((*this).cwiseAbs2().sum()); 00096 } 00097 00104 template<typename Derived> 00105 inline typename NumTraits<typename internal::traits<Derived>::Scalar>::Real MatrixBase<Derived>::norm() const 00106 { 00107 return numext::sqrt(squaredNorm()); 00108 } 00109 00119 template<typename Derived> 00120 inline const typename MatrixBase<Derived>::PlainObject 00121 MatrixBase<Derived>::normalized() const 00122 { 00123 typedef typename internal::nested_eval<Derived,2>::type _Nested; 00124 _Nested n(derived()); 00125 RealScalar z = n.squaredNorm(); 00126 // NOTE: after extensive benchmarking, this conditional does not impact performance, at least on recent x86 CPU 00127 if(z>RealScalar(0)) 00128 return n / numext::sqrt(z); 00129 else 00130 return n; 00131 } 00132 00141 template<typename Derived> 00142 inline void MatrixBase<Derived>::normalize() 00143 { 00144 RealScalar z = squaredNorm(); 00145 // NOTE: after extensive benchmarking, this conditional does not impact performance, at least on recent x86 CPU 00146 if(z>RealScalar(0)) 00147 derived() /= numext::sqrt(z); 00148 } 00149 00162 template<typename Derived> 00163 inline const typename MatrixBase<Derived>::PlainObject 00164 MatrixBase<Derived>::stableNormalized() const 00165 { 00166 typedef typename internal::nested_eval<Derived,3>::type _Nested; 00167 _Nested n(derived()); 00168 RealScalar w = n.cwiseAbs().maxCoeff(); 00169 RealScalar z = (n/w).squaredNorm(); 00170 if(z>RealScalar(0)) 00171 return n / (numext::sqrt(z)*w); 00172 else 00173 return n; 00174 } 00175 00187 template<typename Derived> 00188 inline void MatrixBase<Derived>::stableNormalize() 00189 { 00190 RealScalar w = cwiseAbs().maxCoeff(); 00191 RealScalar z = (derived()/w).squaredNorm(); 00192 if(z>RealScalar(0)) 00193 derived() /= numext::sqrt(z)*w; 00194 } 00195 00196 //---------- implementation of other norms ---------- 00197 00198 namespace internal { 00199 00200 template<typename Derived, int p> 00201 struct lpNorm_selector 00202 { 00203 typedef typename NumTraits<typename traits<Derived>::Scalar>::Real RealScalar; 00204 EIGEN_DEVICE_FUNC 00205 static inline RealScalar run(const MatrixBase<Derived>& m) 00206 { 00207 EIGEN_USING_STD_MATH(pow) 00208 return pow(m.cwiseAbs().array().pow(p).sum(), RealScalar(1)/p); 00209 } 00210 }; 00211 00212 template<typename Derived> 00213 struct lpNorm_selector<Derived, 1> 00214 { 00215 EIGEN_DEVICE_FUNC 00216 static inline typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m) 00217 { 00218 return m.cwiseAbs().sum(); 00219 } 00220 }; 00221 00222 template<typename Derived> 00223 struct lpNorm_selector<Derived, 2> 00224 { 00225 EIGEN_DEVICE_FUNC 00226 static inline typename NumTraits<typename traits<Derived>::Scalar>::Real run(const MatrixBase<Derived>& m) 00227 { 00228 return m.norm(); 00229 } 00230 }; 00231 00232 template<typename Derived> 00233 struct lpNorm_selector<Derived, Infinity> 00234 { 00235 typedef typename NumTraits<typename traits<Derived>::Scalar>::Real RealScalar; 00236 EIGEN_DEVICE_FUNC 00237 static inline RealScalar run(const MatrixBase<Derived>& m) 00238 { 00239 if(Derived::SizeAtCompileTime==0 || (Derived::SizeAtCompileTime==Dynamic && m.size()==0)) 00240 return RealScalar(0); 00241 return m.cwiseAbs().maxCoeff(); 00242 } 00243 }; 00244 00245 } // end namespace internal 00246 00257 template<typename Derived> 00258 template<int p> 00259 #ifndef EIGEN_PARSED_BY_DOXYGEN 00260 inline typename NumTraits<typename internal::traits<Derived>::Scalar>::Real 00261 #else 00262 MatrixBase<Derived>::RealScalar 00263 #endif 00264 MatrixBase<Derived>::lpNorm() const 00265 { 00266 return internal::lpNorm_selector<Derived, p>::run(*this); 00267 } 00268 00269 //---------- implementation of isOrthogonal / isUnitary ---------- 00270 00277 template<typename Derived> 00278 template<typename OtherDerived> 00279 bool MatrixBase<Derived>::isOrthogonal 00280 (const MatrixBase<OtherDerived>& other, const RealScalar& prec) const 00281 { 00282 typename internal::nested_eval<Derived,2>::type nested(derived()); 00283 typename internal::nested_eval<OtherDerived,2>::type otherNested(other.derived()); 00284 return numext::abs2(nested.dot(otherNested)) <= prec * prec * nested.squaredNorm() * otherNested.squaredNorm(); 00285 } 00286 00298 template<typename Derived> 00299 bool MatrixBase<Derived>::isUnitary(const RealScalar& prec) const 00300 { 00301 typename internal::nested_eval<Derived,1>::type self(derived()); 00302 for(Index i = 0; i < cols(); ++i) 00303 { 00304 if(!internal::isApprox(self.col(i).squaredNorm(), static_cast<RealScalar>(1), prec)) 00305 return false; 00306 for(Index j = 0; j < i; ++j) 00307 if(!internal::isMuchSmallerThan(self.col(i).dot(self.col(j)), static_cast<Scalar>(1), prec)) 00308 return false; 00309 } 00310 return true; 00311 } 00312 00313 } // end namespace Eigen 00314 00315 #endif // EIGEN_DOT_H