Eigen  3.2.7
SparseCwiseBinaryOp.h
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 // 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_SPARSE_CWISE_BINARY_OP_H
11 #define EIGEN_SPARSE_CWISE_BINARY_OP_H
12 
13 namespace Eigen {
14 
15 // Here we have to handle 3 cases:
16 // 1 - sparse op dense
17 // 2 - dense op sparse
18 // 3 - sparse op sparse
19 // We also need to implement a 4th iterator for:
20 // 4 - dense op dense
21 // Finally, we also need to distinguish between the product and other operations :
22 // configuration returned mode
23 // 1 - sparse op dense product sparse
24 // generic dense
25 // 2 - dense op sparse product sparse
26 // generic dense
27 // 3 - sparse op sparse product sparse
28 // generic sparse
29 // 4 - dense op dense product dense
30 // generic dense
31 
32 namespace internal {
33 
34 template<> struct promote_storage_type<Dense,Sparse>
35 { typedef Sparse ret; };
36 
37 template<> struct promote_storage_type<Sparse,Dense>
38 { typedef Sparse ret; };
39 
40 template<typename BinaryOp, typename Lhs, typename Rhs, typename Derived,
41  typename _LhsStorageMode = typename traits<Lhs>::StorageKind,
42  typename _RhsStorageMode = typename traits<Rhs>::StorageKind>
43 class sparse_cwise_binary_op_inner_iterator_selector;
44 
45 } // end namespace internal
46 
47 template<typename BinaryOp, typename Lhs, typename Rhs>
48 class CwiseBinaryOpImpl<BinaryOp, Lhs, Rhs, Sparse>
49  : public SparseMatrixBase<CwiseBinaryOp<BinaryOp, Lhs, Rhs> >
50 {
51  public:
52  class InnerIterator;
53  class ReverseInnerIterator;
54  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> Derived;
55  EIGEN_SPARSE_PUBLIC_INTERFACE(Derived)
56  CwiseBinaryOpImpl()
57  {
58  typedef typename internal::traits<Lhs>::StorageKind LhsStorageKind;
59  typedef typename internal::traits<Rhs>::StorageKind RhsStorageKind;
60  EIGEN_STATIC_ASSERT((
61  (!internal::is_same<LhsStorageKind,RhsStorageKind>::value)
62  || ((Lhs::Flags&RowMajorBit) == (Rhs::Flags&RowMajorBit))),
63  THE_STORAGE_ORDER_OF_BOTH_SIDES_MUST_MATCH);
64  }
65 };
66 
67 template<typename BinaryOp, typename Lhs, typename Rhs>
68 class CwiseBinaryOpImpl<BinaryOp,Lhs,Rhs,Sparse>::InnerIterator
69  : public internal::sparse_cwise_binary_op_inner_iterator_selector<BinaryOp,Lhs,Rhs,typename CwiseBinaryOpImpl<BinaryOp,Lhs,Rhs,Sparse>::InnerIterator>
70 {
71  public:
72  typedef typename Lhs::Index Index;
73  typedef internal::sparse_cwise_binary_op_inner_iterator_selector<
74  BinaryOp,Lhs,Rhs, InnerIterator> Base;
75 
76  // NOTE: we have to prefix Index by "typename Lhs::" to avoid an ICE with VC11
77  EIGEN_STRONG_INLINE InnerIterator(const CwiseBinaryOpImpl& binOp, typename Lhs::Index outer)
78  : Base(binOp.derived(),outer)
79  {}
80 };
81 
82 /***************************************************************************
83 * Implementation of inner-iterators
84 ***************************************************************************/
85 
86 // template<typename T> struct internal::func_is_conjunction { enum { ret = false }; };
87 // template<typename T> struct internal::func_is_conjunction<internal::scalar_product_op<T> > { enum { ret = true }; };
88 
89 // TODO generalize the internal::scalar_product_op specialization to all conjunctions if any !
90 
91 namespace internal {
92 
93 // sparse - sparse (generic)
94 template<typename BinaryOp, typename Lhs, typename Rhs, typename Derived>
95 class sparse_cwise_binary_op_inner_iterator_selector<BinaryOp, Lhs, Rhs, Derived, Sparse, Sparse>
96 {
97  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> CwiseBinaryXpr;
98  typedef typename traits<CwiseBinaryXpr>::Scalar Scalar;
99  typedef typename traits<CwiseBinaryXpr>::_LhsNested _LhsNested;
100  typedef typename traits<CwiseBinaryXpr>::_RhsNested _RhsNested;
101  typedef typename _LhsNested::InnerIterator LhsIterator;
102  typedef typename _RhsNested::InnerIterator RhsIterator;
103  typedef typename Lhs::Index Index;
104 
105  public:
106 
107  EIGEN_STRONG_INLINE sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, Index outer)
108  : m_lhsIter(xpr.lhs(),outer), m_rhsIter(xpr.rhs(),outer), m_functor(xpr.functor())
109  {
110  this->operator++();
111  }
112 
113  EIGEN_STRONG_INLINE Derived& operator++()
114  {
115  if (m_lhsIter && m_rhsIter && (m_lhsIter.index() == m_rhsIter.index()))
116  {
117  m_id = m_lhsIter.index();
118  m_value = m_functor(m_lhsIter.value(), m_rhsIter.value());
119  ++m_lhsIter;
120  ++m_rhsIter;
121  }
122  else if (m_lhsIter && (!m_rhsIter || (m_lhsIter.index() < m_rhsIter.index())))
123  {
124  m_id = m_lhsIter.index();
125  m_value = m_functor(m_lhsIter.value(), Scalar(0));
126  ++m_lhsIter;
127  }
128  else if (m_rhsIter && (!m_lhsIter || (m_lhsIter.index() > m_rhsIter.index())))
129  {
130  m_id = m_rhsIter.index();
131  m_value = m_functor(Scalar(0), m_rhsIter.value());
132  ++m_rhsIter;
133  }
134  else
135  {
136  m_value = 0; // this is to avoid a compilation warning
137  m_id = -1;
138  }
139  return *static_cast<Derived*>(this);
140  }
141 
142  EIGEN_STRONG_INLINE Scalar value() const { return m_value; }
143 
144  EIGEN_STRONG_INLINE Index index() const { return m_id; }
145  EIGEN_STRONG_INLINE Index row() const { return Lhs::IsRowMajor ? m_lhsIter.row() : index(); }
146  EIGEN_STRONG_INLINE Index col() const { return Lhs::IsRowMajor ? index() : m_lhsIter.col(); }
147 
148  EIGEN_STRONG_INLINE operator bool() const { return m_id>=0; }
149 
150  protected:
151  LhsIterator m_lhsIter;
152  RhsIterator m_rhsIter;
153  const BinaryOp& m_functor;
154  Scalar m_value;
155  Index m_id;
156 };
157 
158 // sparse - sparse (product)
159 template<typename T, typename Lhs, typename Rhs, typename Derived>
160 class sparse_cwise_binary_op_inner_iterator_selector<scalar_product_op<T>, Lhs, Rhs, Derived, Sparse, Sparse>
161 {
162  typedef scalar_product_op<T> BinaryFunc;
163  typedef CwiseBinaryOp<BinaryFunc, Lhs, Rhs> CwiseBinaryXpr;
164  typedef typename CwiseBinaryXpr::Scalar Scalar;
165  typedef typename traits<CwiseBinaryXpr>::_LhsNested _LhsNested;
166  typedef typename _LhsNested::InnerIterator LhsIterator;
167  typedef typename traits<CwiseBinaryXpr>::_RhsNested _RhsNested;
168  typedef typename _RhsNested::InnerIterator RhsIterator;
169  typedef typename Lhs::Index Index;
170  public:
171 
172  EIGEN_STRONG_INLINE sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, Index outer)
173  : m_lhsIter(xpr.lhs(),outer), m_rhsIter(xpr.rhs(),outer), m_functor(xpr.functor())
174  {
175  while (m_lhsIter && m_rhsIter && (m_lhsIter.index() != m_rhsIter.index()))
176  {
177  if (m_lhsIter.index() < m_rhsIter.index())
178  ++m_lhsIter;
179  else
180  ++m_rhsIter;
181  }
182  }
183 
184  EIGEN_STRONG_INLINE Derived& operator++()
185  {
186  ++m_lhsIter;
187  ++m_rhsIter;
188  while (m_lhsIter && m_rhsIter && (m_lhsIter.index() != m_rhsIter.index()))
189  {
190  if (m_lhsIter.index() < m_rhsIter.index())
191  ++m_lhsIter;
192  else
193  ++m_rhsIter;
194  }
195  return *static_cast<Derived*>(this);
196  }
197 
198  EIGEN_STRONG_INLINE Scalar value() const { return m_functor(m_lhsIter.value(), m_rhsIter.value()); }
199 
200  EIGEN_STRONG_INLINE Index index() const { return m_lhsIter.index(); }
201  EIGEN_STRONG_INLINE Index row() const { return m_lhsIter.row(); }
202  EIGEN_STRONG_INLINE Index col() const { return m_lhsIter.col(); }
203 
204  EIGEN_STRONG_INLINE operator bool() const { return (m_lhsIter && m_rhsIter); }
205 
206  protected:
207  LhsIterator m_lhsIter;
208  RhsIterator m_rhsIter;
209  const BinaryFunc& m_functor;
210 };
211 
212 // sparse - dense (product)
213 template<typename T, typename Lhs, typename Rhs, typename Derived>
214 class sparse_cwise_binary_op_inner_iterator_selector<scalar_product_op<T>, Lhs, Rhs, Derived, Sparse, Dense>
215 {
216  typedef scalar_product_op<T> BinaryFunc;
217  typedef CwiseBinaryOp<BinaryFunc, Lhs, Rhs> CwiseBinaryXpr;
218  typedef typename CwiseBinaryXpr::Scalar Scalar;
219  typedef typename traits<CwiseBinaryXpr>::_LhsNested _LhsNested;
220  typedef typename traits<CwiseBinaryXpr>::RhsNested RhsNested;
221  typedef typename _LhsNested::InnerIterator LhsIterator;
222  typedef typename Lhs::Index Index;
223  enum { IsRowMajor = (int(Lhs::Flags)&RowMajorBit)==RowMajorBit };
224  public:
225 
226  EIGEN_STRONG_INLINE sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, Index outer)
227  : m_rhs(xpr.rhs()), m_lhsIter(xpr.lhs(),outer), m_functor(xpr.functor()), m_outer(outer)
228  {}
229 
230  EIGEN_STRONG_INLINE Derived& operator++()
231  {
232  ++m_lhsIter;
233  return *static_cast<Derived*>(this);
234  }
235 
236  EIGEN_STRONG_INLINE Scalar value() const
237  { return m_functor(m_lhsIter.value(),
238  m_rhs.coeff(IsRowMajor?m_outer:m_lhsIter.index(),IsRowMajor?m_lhsIter.index():m_outer)); }
239 
240  EIGEN_STRONG_INLINE Index index() const { return m_lhsIter.index(); }
241  EIGEN_STRONG_INLINE Index row() const { return m_lhsIter.row(); }
242  EIGEN_STRONG_INLINE Index col() const { return m_lhsIter.col(); }
243 
244  EIGEN_STRONG_INLINE operator bool() const { return m_lhsIter; }
245 
246  protected:
247  RhsNested m_rhs;
248  LhsIterator m_lhsIter;
249  const BinaryFunc m_functor;
250  const Index m_outer;
251 };
252 
253 // sparse - dense (product)
254 template<typename T, typename Lhs, typename Rhs, typename Derived>
255 class sparse_cwise_binary_op_inner_iterator_selector<scalar_product_op<T>, Lhs, Rhs, Derived, Dense, Sparse>
256 {
257  typedef scalar_product_op<T> BinaryFunc;
258  typedef CwiseBinaryOp<BinaryFunc, Lhs, Rhs> CwiseBinaryXpr;
259  typedef typename CwiseBinaryXpr::Scalar Scalar;
260  typedef typename traits<CwiseBinaryXpr>::_RhsNested _RhsNested;
261  typedef typename _RhsNested::InnerIterator RhsIterator;
262  typedef typename Lhs::Index Index;
263 
264  enum { IsRowMajor = (int(Rhs::Flags)&RowMajorBit)==RowMajorBit };
265  public:
266 
267  EIGEN_STRONG_INLINE sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, Index outer)
268  : m_xpr(xpr), m_rhsIter(xpr.rhs(),outer), m_functor(xpr.functor()), m_outer(outer)
269  {}
270 
271  EIGEN_STRONG_INLINE Derived& operator++()
272  {
273  ++m_rhsIter;
274  return *static_cast<Derived*>(this);
275  }
276 
277  EIGEN_STRONG_INLINE Scalar value() const
278  { return m_functor(m_xpr.lhs().coeff(IsRowMajor?m_outer:m_rhsIter.index(),IsRowMajor?m_rhsIter.index():m_outer), m_rhsIter.value()); }
279 
280  EIGEN_STRONG_INLINE Index index() const { return m_rhsIter.index(); }
281  EIGEN_STRONG_INLINE Index row() const { return m_rhsIter.row(); }
282  EIGEN_STRONG_INLINE Index col() const { return m_rhsIter.col(); }
283 
284  EIGEN_STRONG_INLINE operator bool() const { return m_rhsIter; }
285 
286  protected:
287  const CwiseBinaryXpr& m_xpr;
288  RhsIterator m_rhsIter;
289  const BinaryFunc& m_functor;
290  const Index m_outer;
291 };
292 
293 } // end namespace internal
294 
295 /***************************************************************************
296 * Implementation of SparseMatrixBase and SparseCwise functions/operators
297 ***************************************************************************/
298 
299 template<typename Derived>
300 template<typename OtherDerived>
301 EIGEN_STRONG_INLINE Derived &
302 SparseMatrixBase<Derived>::operator-=(const SparseMatrixBase<OtherDerived> &other)
303 {
304  return derived() = derived() - other.derived();
305 }
306 
307 template<typename Derived>
308 template<typename OtherDerived>
309 EIGEN_STRONG_INLINE Derived &
310 SparseMatrixBase<Derived>::operator+=(const SparseMatrixBase<OtherDerived>& other)
311 {
312  return derived() = derived() + other.derived();
313 }
314 
315 template<typename Derived>
316 template<typename OtherDerived>
317 EIGEN_STRONG_INLINE const typename SparseMatrixBase<Derived>::template CwiseProductDenseReturnType<OtherDerived>::Type
318 SparseMatrixBase<Derived>::cwiseProduct(const MatrixBase<OtherDerived> &other) const
319 {
320  return typename CwiseProductDenseReturnType<OtherDerived>::Type(derived(), other.derived());
321 }
322 
323 } // end namespace Eigen
324 
325 #endif // EIGEN_SPARSE_CWISE_BINARY_OP_H
Definition: LDLT.h:16
const CwiseBinaryOp< internal::scalar_product_op< typename Derived::Scalar, typename OtherDerived::Scalar >, const Derived, const OtherDerived > cwiseProduct(const Eigen::SparseMatrixBase< OtherDerived > &other) const
Definition: SparseMatrixBase.h:23
Definition: Eigen_Colamd.h:54
const unsigned int RowMajorBit
Definition: Constants.h:53