Eigen  3.2.91
SparseCwiseBinaryOp.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008-2014 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 template<typename BinaryOp, typename Lhs, typename Rhs>
33 class CwiseBinaryOpImpl<BinaryOp, Lhs, Rhs, Sparse>
34  : public SparseMatrixBase<CwiseBinaryOp<BinaryOp, Lhs, Rhs> >
35 {
36  public:
37  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> Derived;
38  EIGEN_SPARSE_PUBLIC_INTERFACE(Derived)
39  CwiseBinaryOpImpl()
40  {
41  typedef typename internal::traits<Lhs>::StorageKind LhsStorageKind;
42  typedef typename internal::traits<Rhs>::StorageKind RhsStorageKind;
43  EIGEN_STATIC_ASSERT((
44  (!internal::is_same<LhsStorageKind,RhsStorageKind>::value)
45  || ((Lhs::Flags&RowMajorBit) == (Rhs::Flags&RowMajorBit))),
46  THE_STORAGE_ORDER_OF_BOTH_SIDES_MUST_MATCH);
47  }
48 };
49 
50 namespace internal {
51 
52 template<typename BinaryOp, typename Lhs, typename Rhs, typename Derived,
53  typename _LhsStorageMode = typename traits<Lhs>::StorageKind,
54  typename _RhsStorageMode = typename traits<Rhs>::StorageKind>
55 class sparse_cwise_binary_op_inner_iterator_selector;
56 
57 } // end namespace internal
58 
59 namespace internal {
60 
61 
62 // Generic "sparse OP sparse"
63 template<typename BinaryOp, typename Lhs, typename Rhs>
64 struct binary_evaluator<CwiseBinaryOp<BinaryOp, Lhs, Rhs>, IteratorBased, IteratorBased>
65  : evaluator_base<CwiseBinaryOp<BinaryOp, Lhs, Rhs> >
66 {
67 protected:
68  typedef typename evaluator<Lhs>::InnerIterator LhsIterator;
69  typedef typename evaluator<Rhs>::InnerIterator RhsIterator;
70  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> XprType;
71  typedef typename traits<XprType>::Scalar Scalar;
72  typedef typename XprType::StorageIndex StorageIndex;
73 public:
74 
75  class ReverseInnerIterator;
76  class InnerIterator
77  {
78  public:
79 
80  EIGEN_STRONG_INLINE InnerIterator(const binary_evaluator& aEval, Index outer)
81  : m_lhsIter(aEval.m_lhsImpl,outer), m_rhsIter(aEval.m_rhsImpl,outer), m_functor(aEval.m_functor)
82  {
83  this->operator++();
84  }
85 
86  EIGEN_STRONG_INLINE InnerIterator& operator++()
87  {
88  if (m_lhsIter && m_rhsIter && (m_lhsIter.index() == m_rhsIter.index()))
89  {
90  m_id = m_lhsIter.index();
91  m_value = m_functor(m_lhsIter.value(), m_rhsIter.value());
92  ++m_lhsIter;
93  ++m_rhsIter;
94  }
95  else if (m_lhsIter && (!m_rhsIter || (m_lhsIter.index() < m_rhsIter.index())))
96  {
97  m_id = m_lhsIter.index();
98  m_value = m_functor(m_lhsIter.value(), Scalar(0));
99  ++m_lhsIter;
100  }
101  else if (m_rhsIter && (!m_lhsIter || (m_lhsIter.index() > m_rhsIter.index())))
102  {
103  m_id = m_rhsIter.index();
104  m_value = m_functor(Scalar(0), m_rhsIter.value());
105  ++m_rhsIter;
106  }
107  else
108  {
109  m_value = 0; // this is to avoid a compilation warning
110  m_id = -1;
111  }
112  return *this;
113  }
114 
115  EIGEN_STRONG_INLINE Scalar value() const { return m_value; }
116 
117  EIGEN_STRONG_INLINE StorageIndex index() const { return m_id; }
118  EIGEN_STRONG_INLINE Index row() const { return Lhs::IsRowMajor ? m_lhsIter.row() : index(); }
119  EIGEN_STRONG_INLINE Index col() const { return Lhs::IsRowMajor ? index() : m_lhsIter.col(); }
120 
121  EIGEN_STRONG_INLINE operator bool() const { return m_id>=0; }
122 
123  protected:
124  LhsIterator m_lhsIter;
125  RhsIterator m_rhsIter;
126  const BinaryOp& m_functor;
127  Scalar m_value;
128  StorageIndex m_id;
129  };
130 
131 
132  enum {
133  CoeffReadCost = evaluator<Lhs>::CoeffReadCost + evaluator<Rhs>::CoeffReadCost + functor_traits<BinaryOp>::Cost,
134  Flags = XprType::Flags
135  };
136 
137  explicit binary_evaluator(const XprType& xpr)
138  : m_functor(xpr.functor()),
139  m_lhsImpl(xpr.lhs()),
140  m_rhsImpl(xpr.rhs())
141  { }
142 
143  inline Index nonZerosEstimate() const {
144  return m_lhsImpl.nonZerosEstimate() + m_rhsImpl.nonZerosEstimate();
145  }
146 
147 protected:
148  const BinaryOp m_functor;
149  evaluator<Lhs> m_lhsImpl;
150  evaluator<Rhs> m_rhsImpl;
151 };
152 
153 // "sparse .* sparse"
154 template<typename T, typename Lhs, typename Rhs>
155 struct binary_evaluator<CwiseBinaryOp<scalar_product_op<T>, Lhs, Rhs>, IteratorBased, IteratorBased>
156  : evaluator_base<CwiseBinaryOp<scalar_product_op<T>, Lhs, Rhs> >
157 {
158 protected:
159  typedef scalar_product_op<T> BinaryOp;
160  typedef typename evaluator<Lhs>::InnerIterator LhsIterator;
161  typedef typename evaluator<Rhs>::InnerIterator RhsIterator;
162  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> XprType;
163  typedef typename XprType::StorageIndex StorageIndex;
164  typedef typename traits<XprType>::Scalar Scalar;
165 public:
166 
167  class ReverseInnerIterator;
168  class InnerIterator
169  {
170  public:
171 
172  EIGEN_STRONG_INLINE InnerIterator(const binary_evaluator& aEval, Index outer)
173  : m_lhsIter(aEval.m_lhsImpl,outer), m_rhsIter(aEval.m_rhsImpl,outer), m_functor(aEval.m_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 InnerIterator& 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 *this;
196  }
197 
198  EIGEN_STRONG_INLINE Scalar value() const { return m_functor(m_lhsIter.value(), m_rhsIter.value()); }
199 
200  EIGEN_STRONG_INLINE StorageIndex 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 BinaryOp& m_functor;
210  };
211 
212 
213  enum {
214  CoeffReadCost = evaluator<Lhs>::CoeffReadCost + evaluator<Rhs>::CoeffReadCost + functor_traits<BinaryOp>::Cost,
215  Flags = XprType::Flags
216  };
217 
218  explicit binary_evaluator(const XprType& xpr)
219  : m_functor(xpr.functor()),
220  m_lhsImpl(xpr.lhs()),
221  m_rhsImpl(xpr.rhs())
222  { }
223 
224  inline Index nonZerosEstimate() const {
225  return (std::min)(m_lhsImpl.nonZerosEstimate(), m_rhsImpl.nonZerosEstimate());
226  }
227 
228 protected:
229  const BinaryOp m_functor;
230  evaluator<Lhs> m_lhsImpl;
231  evaluator<Rhs> m_rhsImpl;
232 };
233 
234 // "dense .* sparse"
235 template<typename T, typename Lhs, typename Rhs>
236 struct binary_evaluator<CwiseBinaryOp<scalar_product_op<T>, Lhs, Rhs>, IndexBased, IteratorBased>
237  : evaluator_base<CwiseBinaryOp<scalar_product_op<T>, Lhs, Rhs> >
238 {
239 protected:
240  typedef scalar_product_op<T> BinaryOp;
241  typedef evaluator<Lhs> LhsEvaluator;
242  typedef typename evaluator<Rhs>::InnerIterator RhsIterator;
243  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> XprType;
244  typedef typename XprType::StorageIndex StorageIndex;
245  typedef typename traits<XprType>::Scalar Scalar;
246 public:
247 
248  class ReverseInnerIterator;
249  class InnerIterator
250  {
251  enum { IsRowMajor = (int(Rhs::Flags)&RowMajorBit)==RowMajorBit };
252 
253  public:
254 
255  EIGEN_STRONG_INLINE InnerIterator(const binary_evaluator& aEval, Index outer)
256  : m_lhsEval(aEval.m_lhsImpl), m_rhsIter(aEval.m_rhsImpl,outer), m_functor(aEval.m_functor), m_outer(outer)
257  {}
258 
259  EIGEN_STRONG_INLINE InnerIterator& operator++()
260  {
261  ++m_rhsIter;
262  return *this;
263  }
264 
265  EIGEN_STRONG_INLINE Scalar value() const
266  { return m_functor(m_lhsEval.coeff(IsRowMajor?m_outer:m_rhsIter.index(),IsRowMajor?m_rhsIter.index():m_outer), m_rhsIter.value()); }
267 
268  EIGEN_STRONG_INLINE StorageIndex index() const { return m_rhsIter.index(); }
269  EIGEN_STRONG_INLINE Index row() const { return m_rhsIter.row(); }
270  EIGEN_STRONG_INLINE Index col() const { return m_rhsIter.col(); }
271 
272  EIGEN_STRONG_INLINE operator bool() const { return m_rhsIter; }
273 
274  protected:
275  const LhsEvaluator &m_lhsEval;
276  RhsIterator m_rhsIter;
277  const BinaryOp& m_functor;
278  const Index m_outer;
279  };
280 
281 
282  enum {
283  CoeffReadCost = evaluator<Lhs>::CoeffReadCost + evaluator<Rhs>::CoeffReadCost + functor_traits<BinaryOp>::Cost,
284  Flags = XprType::Flags
285  };
286 
287  explicit binary_evaluator(const XprType& xpr)
288  : m_functor(xpr.functor()),
289  m_lhsImpl(xpr.lhs()),
290  m_rhsImpl(xpr.rhs())
291  { }
292 
293  inline Index nonZerosEstimate() const {
294  return m_rhsImpl.nonZerosEstimate();
295  }
296 
297 protected:
298  const BinaryOp m_functor;
299  evaluator<Lhs> m_lhsImpl;
300  evaluator<Rhs> m_rhsImpl;
301 };
302 
303 // "sparse .* dense"
304 template<typename T, typename Lhs, typename Rhs>
305 struct binary_evaluator<CwiseBinaryOp<scalar_product_op<T>, Lhs, Rhs>, IteratorBased, IndexBased>
306  : evaluator_base<CwiseBinaryOp<scalar_product_op<T>, Lhs, Rhs> >
307 {
308 protected:
309  typedef scalar_product_op<T> BinaryOp;
310  typedef typename evaluator<Lhs>::InnerIterator LhsIterator;
311  typedef evaluator<Rhs> RhsEvaluator;
312  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> XprType;
313  typedef typename XprType::StorageIndex StorageIndex;
314  typedef typename traits<XprType>::Scalar Scalar;
315 public:
316 
317  class ReverseInnerIterator;
318  class InnerIterator
319  {
320  enum { IsRowMajor = (int(Lhs::Flags)&RowMajorBit)==RowMajorBit };
321 
322  public:
323 
324  EIGEN_STRONG_INLINE InnerIterator(const binary_evaluator& aEval, Index outer)
325  : m_lhsIter(aEval.m_lhsImpl,outer), m_rhsEval(aEval.m_rhsImpl), m_functor(aEval.m_functor), m_outer(outer)
326  {}
327 
328  EIGEN_STRONG_INLINE InnerIterator& operator++()
329  {
330  ++m_lhsIter;
331  return *this;
332  }
333 
334  EIGEN_STRONG_INLINE Scalar value() const
335  { return m_functor(m_lhsIter.value(),
336  m_rhsEval.coeff(IsRowMajor?m_outer:m_lhsIter.index(),IsRowMajor?m_lhsIter.index():m_outer)); }
337 
338  EIGEN_STRONG_INLINE StorageIndex index() const { return m_lhsIter.index(); }
339  EIGEN_STRONG_INLINE Index row() const { return m_lhsIter.row(); }
340  EIGEN_STRONG_INLINE Index col() const { return m_lhsIter.col(); }
341 
342  EIGEN_STRONG_INLINE operator bool() const { return m_lhsIter; }
343 
344  protected:
345  LhsIterator m_lhsIter;
346  const evaluator<Rhs> &m_rhsEval;
347  const BinaryOp& m_functor;
348  const Index m_outer;
349  };
350 
351 
352  enum {
353  CoeffReadCost = evaluator<Lhs>::CoeffReadCost + evaluator<Rhs>::CoeffReadCost + functor_traits<BinaryOp>::Cost,
354  Flags = XprType::Flags
355  };
356 
357  explicit binary_evaluator(const XprType& xpr)
358  : m_functor(xpr.functor()),
359  m_lhsImpl(xpr.lhs()),
360  m_rhsImpl(xpr.rhs())
361  { }
362 
363  inline Index nonZerosEstimate() const {
364  return m_lhsImpl.nonZerosEstimate();
365  }
366 
367 protected:
368  const BinaryOp m_functor;
369  evaluator<Lhs> m_lhsImpl;
370  evaluator<Rhs> m_rhsImpl;
371 };
372 
373 }
374 
375 /***************************************************************************
376 * Implementation of SparseMatrixBase and SparseCwise functions/operators
377 ***************************************************************************/
378 
379 template<typename Derived>
380 template<typename OtherDerived>
381 EIGEN_STRONG_INLINE Derived &
382 SparseMatrixBase<Derived>::operator-=(const SparseMatrixBase<OtherDerived> &other)
383 {
384  return derived() = derived() - other.derived();
385 }
386 
387 template<typename Derived>
388 template<typename OtherDerived>
389 EIGEN_STRONG_INLINE Derived &
390 SparseMatrixBase<Derived>::operator+=(const SparseMatrixBase<OtherDerived>& other)
391 {
392  return derived() = derived() + other.derived();
393 }
394 
395 template<typename Derived>
396 template<typename OtherDerived>
397 Derived& SparseMatrixBase<Derived>::operator+=(const DiagonalBase<OtherDerived>& other)
398 {
399  call_assignment_no_alias(derived(), other.derived(), internal::add_assign_op<Scalar>());
400  return derived();
401 }
402 
403 template<typename Derived>
404 template<typename OtherDerived>
405 Derived& SparseMatrixBase<Derived>::operator-=(const DiagonalBase<OtherDerived>& other)
406 {
407  call_assignment_no_alias(derived(), other.derived(), internal::sub_assign_op<Scalar>());
408  return derived();
409 }
410 
411 template<typename Derived>
412 template<typename OtherDerived>
413 EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE
414 SparseMatrixBase<Derived>::cwiseProduct(const MatrixBase<OtherDerived> &other) const
415 {
416  return EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE(derived(), other.derived());
417 }
418 
419 } // end namespace Eigen
420 
421 #endif // EIGEN_SPARSE_CWISE_BINARY_OP_H
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:24
Definition: LDLT.h:16
const unsigned int RowMajorBit
Definition: Constants.h:53
Definition: Eigen_Colamd.h:54