SparseDot.h
Go to the documentation of this file.
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_SPARSE_DOT_H
26 #define EIGEN_SPARSE_DOT_H
27 
28 namespace Eigen {
29 
30 template<typename Derived>
31 template<typename OtherDerived>
32 typename internal::traits<Derived>::Scalar
34 {
37  EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
38  EIGEN_STATIC_ASSERT((internal::is_same<Scalar, typename OtherDerived::Scalar>::value),
39  YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
40 
41  eigen_assert(size() == other.size());
42  eigen_assert(other.size()>0 && "you are using a non initialized vector");
43 
44  typename Derived::InnerIterator i(derived(),0);
45  Scalar res(0);
46  while (i)
47  {
48  res += internal::conj(i.value()) * other.coeff(i.index());
49  ++i;
50  }
51  return res;
52 }
53 
54 template<typename Derived>
55 template<typename OtherDerived>
56 typename internal::traits<Derived>::Scalar
58 {
61  EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
62  EIGEN_STATIC_ASSERT((internal::is_same<Scalar, typename OtherDerived::Scalar>::value),
63  YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
64 
65  eigen_assert(size() == other.size());
66 
67  typedef typename Derived::Nested Nested;
68  typedef typename OtherDerived::Nested OtherNested;
69  typedef typename internal::remove_all<Nested>::type NestedCleaned;
70  typedef typename internal::remove_all<OtherNested>::type OtherNestedCleaned;
71 
72  const Nested nthis(derived());
73  const OtherNested nother(other.derived());
74 
75  typename NestedCleaned::InnerIterator i(nthis,0);
76  typename OtherNestedCleaned::InnerIterator j(nother,0);
77  Scalar res(0);
78  while (i && j)
79  {
80  if (i.index()==j.index())
81  {
82  res += internal::conj(i.value()) * j.value();
83  ++i; ++j;
84  }
85  else if (i.index()<j.index())
86  ++i;
87  else
88  ++j;
89  }
90  return res;
91 }
92 
93 template<typename Derived>
96 {
97  return internal::real((*this).cwiseAbs2().sum());
98 }
99 
100 template<typename Derived>
103 {
104  return internal::sqrt(squaredNorm());
105 }
106 
107 } // end namespace Eigen
108 
109 #endif // EIGEN_SPARSE_DOT_H