$darkmode
Eigen  5.0.1-dev
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 // IWYU pragma: private
14 #include "./InternalHeaderCheck.h"
15 
16 namespace Eigen {
17 
18 // Here we have to handle 3 cases:
19 // 1 - sparse op dense
20 // 2 - dense op sparse
21 // 3 - sparse op sparse
22 // We also need to implement a 4th iterator for:
23 // 4 - dense op dense
24 // Finally, we also need to distinguish between the product and other operations :
25 // configuration returned mode
26 // 1 - sparse op dense product sparse
27 // generic dense
28 // 2 - dense op sparse product sparse
29 // generic dense
30 // 3 - sparse op sparse product sparse
31 // generic sparse
32 // 4 - dense op dense product dense
33 // generic dense
34 //
35 // TODO to ease compiler job, we could specialize product/quotient with a scalar
36 // and fallback to cwise-unary evaluator using bind1st_op and bind2nd_op.
37 
38 template <typename BinaryOp, typename Lhs, typename Rhs>
39 class CwiseBinaryOpImpl<BinaryOp, Lhs, Rhs, Sparse> : public SparseMatrixBase<CwiseBinaryOp<BinaryOp, Lhs, Rhs> > {
40  public:
41  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> Derived;
42  typedef SparseMatrixBase<Derived> Base;
43  EIGEN_SPARSE_PUBLIC_INTERFACE(Derived)
44  EIGEN_STATIC_ASSERT(((!internal::is_same<typename internal::traits<Lhs>::StorageKind,
45  typename internal::traits<Rhs>::StorageKind>::value) ||
46  ((internal::evaluator<Lhs>::Flags & RowMajorBit) ==
47  (internal::evaluator<Rhs>::Flags & RowMajorBit))),
48  THE_STORAGE_ORDER_OF_BOTH_SIDES_MUST_MATCH)
49 };
50 
51 namespace internal {
52 
53 // The default evaluator performs an "arithmetic" operation on two input arrays.
54 // Given input arrays 'lhs' and 'rhs' and binary functor 'func',
55 // the sparse destination array 'dst' is evaluated as follows:
56 // if lhs(i,j) and rhs(i,j) are present, dst(i,j) = func(lhs(i,j), rhs(i,j))
57 // if lhs(i,j) is present and rhs(i,j) is null, dst(i,j) = func(lhs(i,j), 0)
58 // if lhs(i,j) is null and rhs(i,j) is present, dst(i,j) = func(0, rhs(i,j))
59 
60 // Generic "sparse OP sparse"
61 template <typename XprType>
62 struct binary_sparse_evaluator;
63 
64 template <typename BinaryOp, typename Lhs, typename Rhs>
65 struct binary_evaluator<CwiseBinaryOp<BinaryOp, Lhs, Rhs>, IteratorBased, IteratorBased>
66  : evaluator_base<CwiseBinaryOp<BinaryOp, Lhs, Rhs> > {
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 
74  public:
75  class InnerIterator {
76  public:
77  EIGEN_STRONG_INLINE InnerIterator(const binary_evaluator& aEval, Index outer)
78  : m_lhsIter(aEval.m_lhsImpl, outer),
79  m_rhsIter(aEval.m_rhsImpl, outer),
80  m_functor(aEval.m_functor),
81  m_value(Scalar(0)) {
82  this->operator++();
83  }
84 
85  EIGEN_STRONG_INLINE InnerIterator& operator++() {
86  if (m_lhsIter && m_rhsIter && (m_lhsIter.index() == m_rhsIter.index())) {
87  m_id = m_lhsIter.index();
88  m_value = m_functor(m_lhsIter.value(), m_rhsIter.value());
89  ++m_lhsIter;
90  ++m_rhsIter;
91  } else if (m_lhsIter && (!m_rhsIter || (m_lhsIter.index() < m_rhsIter.index()))) {
92  m_id = m_lhsIter.index();
93  m_value = m_functor(m_lhsIter.value(), Scalar(0));
94  ++m_lhsIter;
95  } else if (m_rhsIter && (!m_lhsIter || (m_lhsIter.index() > m_rhsIter.index()))) {
96  m_id = m_rhsIter.index();
97  m_value = m_functor(Scalar(0), m_rhsIter.value());
98  ++m_rhsIter;
99  } else {
100  m_id = -1;
101  }
102  return *this;
103  }
104 
105  EIGEN_STRONG_INLINE Scalar value() const { return m_value; }
106 
107  EIGEN_STRONG_INLINE StorageIndex index() const { return m_id; }
108  EIGEN_STRONG_INLINE Index outer() const { return m_lhsIter.outer(); }
109  EIGEN_STRONG_INLINE Index row() const { return Lhs::IsRowMajor ? m_lhsIter.row() : index(); }
110  EIGEN_STRONG_INLINE Index col() const { return Lhs::IsRowMajor ? index() : m_lhsIter.col(); }
111 
112  EIGEN_STRONG_INLINE operator bool() const { return m_id >= 0; }
113 
114  protected:
115  LhsIterator m_lhsIter;
116  RhsIterator m_rhsIter;
117  const BinaryOp& m_functor;
118  Scalar m_value;
119  StorageIndex m_id;
120  };
121 
122  enum {
123  CoeffReadCost =
124  int(evaluator<Lhs>::CoeffReadCost) + int(evaluator<Rhs>::CoeffReadCost) + int(functor_traits<BinaryOp>::Cost),
125  Flags = XprType::Flags
126  };
127 
128  explicit binary_evaluator(const XprType& xpr) : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()) {
129  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
130  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
131  }
132 
133  inline Index nonZerosEstimate() const { return m_lhsImpl.nonZerosEstimate() + m_rhsImpl.nonZerosEstimate(); }
134 
135  protected:
136  const BinaryOp m_functor;
137  evaluator<Lhs> m_lhsImpl;
138  evaluator<Rhs> m_rhsImpl;
139 };
140 
141 // dense op sparse
142 template <typename BinaryOp, typename Lhs, typename Rhs>
143 struct binary_evaluator<CwiseBinaryOp<BinaryOp, Lhs, Rhs>, IndexBased, IteratorBased>
144  : evaluator_base<CwiseBinaryOp<BinaryOp, Lhs, Rhs> > {
145  protected:
146  typedef typename evaluator<Rhs>::InnerIterator RhsIterator;
147  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> XprType;
148  typedef typename traits<XprType>::Scalar Scalar;
149  typedef typename XprType::StorageIndex StorageIndex;
150 
151  public:
152  class InnerIterator {
153  enum { IsRowMajor = (int(Rhs::Flags) & RowMajorBit) == RowMajorBit };
154 
155  public:
156  EIGEN_STRONG_INLINE InnerIterator(const binary_evaluator& aEval, Index outer)
157  : m_lhsEval(aEval.m_lhsImpl),
158  m_rhsIter(aEval.m_rhsImpl, outer),
159  m_functor(aEval.m_functor),
160  m_value(0),
161  m_id(-1),
162  m_innerSize(aEval.m_expr.rhs().innerSize()) {
163  this->operator++();
164  }
165 
166  EIGEN_STRONG_INLINE InnerIterator& operator++() {
167  ++m_id;
168  if (m_id < m_innerSize) {
169  Scalar lhsVal = m_lhsEval.coeff(IsRowMajor ? m_rhsIter.outer() : m_id, IsRowMajor ? m_id : m_rhsIter.outer());
170  if (m_rhsIter && m_rhsIter.index() == m_id) {
171  m_value = m_functor(lhsVal, m_rhsIter.value());
172  ++m_rhsIter;
173  } else
174  m_value = m_functor(lhsVal, Scalar(0));
175  }
176 
177  return *this;
178  }
179 
180  EIGEN_STRONG_INLINE Scalar value() const {
181  eigen_internal_assert(m_id < m_innerSize);
182  return m_value;
183  }
184 
185  EIGEN_STRONG_INLINE StorageIndex index() const { return m_id; }
186  EIGEN_STRONG_INLINE Index outer() const { return m_rhsIter.outer(); }
187  EIGEN_STRONG_INLINE Index row() const { return IsRowMajor ? m_rhsIter.outer() : m_id; }
188  EIGEN_STRONG_INLINE Index col() const { return IsRowMajor ? m_id : m_rhsIter.outer(); }
189 
190  EIGEN_STRONG_INLINE operator bool() const { return m_id < m_innerSize; }
191 
192  protected:
193  const evaluator<Lhs>& m_lhsEval;
194  RhsIterator m_rhsIter;
195  const BinaryOp& m_functor;
196  Scalar m_value;
197  StorageIndex m_id;
198  StorageIndex m_innerSize;
199  };
200 
201  enum {
202  CoeffReadCost =
203  int(evaluator<Lhs>::CoeffReadCost) + int(evaluator<Rhs>::CoeffReadCost) + int(functor_traits<BinaryOp>::Cost),
204  Flags = XprType::Flags
205  };
206 
207  explicit binary_evaluator(const XprType& xpr)
208  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()), m_expr(xpr) {
209  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
210  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
211  }
212 
213  inline Index nonZerosEstimate() const { return m_expr.size(); }
214 
215  protected:
216  const BinaryOp m_functor;
217  evaluator<Lhs> m_lhsImpl;
218  evaluator<Rhs> m_rhsImpl;
219  const XprType& m_expr;
220 };
221 
222 // sparse op dense
223 template <typename BinaryOp, typename Lhs, typename Rhs>
224 struct binary_evaluator<CwiseBinaryOp<BinaryOp, Lhs, Rhs>, IteratorBased, IndexBased>
225  : evaluator_base<CwiseBinaryOp<BinaryOp, Lhs, Rhs> > {
226  protected:
227  typedef typename evaluator<Lhs>::InnerIterator LhsIterator;
228  typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> XprType;
229  typedef typename traits<XprType>::Scalar Scalar;
230  typedef typename XprType::StorageIndex StorageIndex;
231 
232  public:
233  class InnerIterator {
234  enum { IsRowMajor = (int(Lhs::Flags) & RowMajorBit) == RowMajorBit };
235 
236  public:
237  EIGEN_STRONG_INLINE InnerIterator(const binary_evaluator& aEval, Index outer)
238  : m_lhsIter(aEval.m_lhsImpl, outer),
239  m_rhsEval(aEval.m_rhsImpl),
240  m_functor(aEval.m_functor),
241  m_value(0),
242  m_id(-1),
243  m_innerSize(aEval.m_expr.lhs().innerSize()) {
244  this->operator++();
245  }
246 
247  EIGEN_STRONG_INLINE InnerIterator& operator++() {
248  ++m_id;
249  if (m_id < m_innerSize) {
250  Scalar rhsVal = m_rhsEval.coeff(IsRowMajor ? m_lhsIter.outer() : m_id, IsRowMajor ? m_id : m_lhsIter.outer());
251  if (m_lhsIter && m_lhsIter.index() == m_id) {
252  m_value = m_functor(m_lhsIter.value(), rhsVal);
253  ++m_lhsIter;
254  } else
255  m_value = m_functor(Scalar(0), rhsVal);
256  }
257 
258  return *this;
259  }
260 
261  EIGEN_STRONG_INLINE Scalar value() const {
262  eigen_internal_assert(m_id < m_innerSize);
263  return m_value;
264  }
265 
266  EIGEN_STRONG_INLINE StorageIndex index() const { return m_id; }
267  EIGEN_STRONG_INLINE Index outer() const { return m_lhsIter.outer(); }
268  EIGEN_STRONG_INLINE Index row() const { return IsRowMajor ? m_lhsIter.outer() : m_id; }
269  EIGEN_STRONG_INLINE Index col() const { return IsRowMajor ? m_id : m_lhsIter.outer(); }
270 
271  EIGEN_STRONG_INLINE operator bool() const { return m_id < m_innerSize; }
272 
273  protected:
274  LhsIterator m_lhsIter;
275  const evaluator<Rhs>& m_rhsEval;
276  const BinaryOp& m_functor;
277  Scalar m_value;
278  StorageIndex m_id;
279  StorageIndex m_innerSize;
280  };
281 
282  enum {
283  CoeffReadCost =
284  int(evaluator<Lhs>::CoeffReadCost) + int(evaluator<Rhs>::CoeffReadCost) + int(functor_traits<BinaryOp>::Cost),
285  Flags = XprType::Flags
286  };
287 
288  explicit binary_evaluator(const XprType& xpr)
289  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()), m_expr(xpr) {
290  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
291  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
292  }
293 
294  inline Index nonZerosEstimate() const { return m_expr.size(); }
295 
296  protected:
297  const BinaryOp m_functor;
298  evaluator<Lhs> m_lhsImpl;
299  evaluator<Rhs> m_rhsImpl;
300  const XprType& m_expr;
301 };
302 
303 template <typename T, typename LhsKind = typename evaluator_traits<typename T::Lhs>::Kind,
304  typename RhsKind = typename evaluator_traits<typename T::Rhs>::Kind,
305  typename LhsScalar = typename traits<typename T::Lhs>::Scalar,
306  typename RhsScalar = typename traits<typename T::Rhs>::Scalar>
307 struct sparse_conjunction_evaluator;
308 
309 // "sparse .* sparse"
310 template <typename T1, typename T2, typename Lhs, typename Rhs>
311 struct binary_evaluator<CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs>, IteratorBased, IteratorBased>
312  : sparse_conjunction_evaluator<CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs> > {
313  typedef CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs> XprType;
314  typedef sparse_conjunction_evaluator<XprType> Base;
315  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
316 };
317 // "dense .* sparse"
318 template <typename T1, typename T2, typename Lhs, typename Rhs>
319 struct binary_evaluator<CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs>, IndexBased, IteratorBased>
320  : sparse_conjunction_evaluator<CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs> > {
321  typedef CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs> XprType;
322  typedef sparse_conjunction_evaluator<XprType> Base;
323  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
324 };
325 // "sparse .* dense"
326 template <typename T1, typename T2, typename Lhs, typename Rhs>
327 struct binary_evaluator<CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs>, IteratorBased, IndexBased>
328  : sparse_conjunction_evaluator<CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs> > {
329  typedef CwiseBinaryOp<scalar_product_op<T1, T2>, Lhs, Rhs> XprType;
330  typedef sparse_conjunction_evaluator<XprType> Base;
331  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
332 };
333 
334 // "sparse ./ dense"
335 template <typename T1, typename T2, typename Lhs, typename Rhs>
336 struct binary_evaluator<CwiseBinaryOp<scalar_quotient_op<T1, T2>, Lhs, Rhs>, IteratorBased, IndexBased>
337  : sparse_conjunction_evaluator<CwiseBinaryOp<scalar_quotient_op<T1, T2>, Lhs, Rhs> > {
338  typedef CwiseBinaryOp<scalar_quotient_op<T1, T2>, Lhs, Rhs> XprType;
339  typedef sparse_conjunction_evaluator<XprType> Base;
340  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
341 };
342 
343 // "sparse && sparse"
344 template <typename Lhs, typename Rhs>
345 struct binary_evaluator<CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs>, IteratorBased, IteratorBased>
346  : sparse_conjunction_evaluator<CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs> > {
347  typedef CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs> XprType;
348  typedef sparse_conjunction_evaluator<XprType> Base;
349  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
350 };
351 // "dense && sparse"
352 template <typename Lhs, typename Rhs>
353 struct binary_evaluator<CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs>, IndexBased, IteratorBased>
354  : sparse_conjunction_evaluator<CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs> > {
355  typedef CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs> XprType;
356  typedef sparse_conjunction_evaluator<XprType> Base;
357  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
358 };
359 // "sparse && dense"
360 template <typename Lhs, typename Rhs>
361 struct binary_evaluator<CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs>, IteratorBased, IndexBased>
362  : sparse_conjunction_evaluator<CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs> > {
363  typedef CwiseBinaryOp<scalar_boolean_and_op<bool>, Lhs, Rhs> XprType;
364  typedef sparse_conjunction_evaluator<XprType> Base;
365  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
366 };
367 
368 // The conjunction "^" evaluator performs a logical "and" or set "intersection" operation on two input arrays.
369 // Given input arrays 'lhs' and 'rhs' and binary functor 'func',
370 // the sparse destination array 'dst' is evaluated as follows:
371 // if lhs(i,j) and rhs(i,j) are present, dst(i,j) = func(lhs(i,j), rhs(i,j))
372 // if lhs(i,j) is present and rhs(i,j) is null, dst(i,j) is null
373 // if lhs(i,j) is null and rhs(i,j) is present, dst(i,j) is null
374 
375 // "sparse ^ sparse"
376 template <typename XprType>
377 struct sparse_conjunction_evaluator<XprType, IteratorBased, IteratorBased> : evaluator_base<XprType> {
378  protected:
379  typedef typename XprType::Functor BinaryOp;
380  typedef typename XprType::Lhs LhsArg;
381  typedef typename XprType::Rhs RhsArg;
382  typedef typename evaluator<LhsArg>::InnerIterator LhsIterator;
383  typedef typename evaluator<RhsArg>::InnerIterator RhsIterator;
384  typedef typename XprType::StorageIndex StorageIndex;
385  typedef typename traits<XprType>::Scalar Scalar;
386 
387  public:
388  class InnerIterator {
389  public:
390  EIGEN_STRONG_INLINE InnerIterator(const sparse_conjunction_evaluator& aEval, Index outer)
391  : m_lhsIter(aEval.m_lhsImpl, outer), m_rhsIter(aEval.m_rhsImpl, outer), m_functor(aEval.m_functor) {
392  while (m_lhsIter && m_rhsIter && (m_lhsIter.index() != m_rhsIter.index())) {
393  if (m_lhsIter.index() < m_rhsIter.index())
394  ++m_lhsIter;
395  else
396  ++m_rhsIter;
397  }
398  }
399 
400  EIGEN_STRONG_INLINE InnerIterator& operator++() {
401  ++m_lhsIter;
402  ++m_rhsIter;
403  while (m_lhsIter && m_rhsIter && (m_lhsIter.index() != m_rhsIter.index())) {
404  if (m_lhsIter.index() < m_rhsIter.index())
405  ++m_lhsIter;
406  else
407  ++m_rhsIter;
408  }
409  return *this;
410  }
411 
412  EIGEN_STRONG_INLINE Scalar value() const { return m_functor(m_lhsIter.value(), m_rhsIter.value()); }
413 
414  EIGEN_STRONG_INLINE StorageIndex index() const { return m_lhsIter.index(); }
415  EIGEN_STRONG_INLINE Index outer() const { return m_lhsIter.outer(); }
416  EIGEN_STRONG_INLINE Index row() const { return m_lhsIter.row(); }
417  EIGEN_STRONG_INLINE Index col() const { return m_lhsIter.col(); }
418 
419  EIGEN_STRONG_INLINE operator bool() const { return (m_lhsIter && m_rhsIter); }
420 
421  protected:
422  LhsIterator m_lhsIter;
423  RhsIterator m_rhsIter;
424  const BinaryOp& m_functor;
425  };
426 
427  enum {
428  CoeffReadCost = int(evaluator<LhsArg>::CoeffReadCost) + int(evaluator<RhsArg>::CoeffReadCost) +
429  int(functor_traits<BinaryOp>::Cost),
430  Flags = XprType::Flags
431  };
432 
433  explicit sparse_conjunction_evaluator(const XprType& xpr)
434  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()) {
435  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
436  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
437  }
438 
439  inline Index nonZerosEstimate() const {
440  return (std::min)(m_lhsImpl.nonZerosEstimate(), m_rhsImpl.nonZerosEstimate());
441  }
442 
443  protected:
444  const BinaryOp m_functor;
445  evaluator<LhsArg> m_lhsImpl;
446  evaluator<RhsArg> m_rhsImpl;
447 };
448 
449 // "dense ^ sparse"
450 template <typename XprType>
451 struct sparse_conjunction_evaluator<XprType, IndexBased, IteratorBased> : evaluator_base<XprType> {
452  protected:
453  typedef typename XprType::Functor BinaryOp;
454  typedef typename XprType::Lhs LhsArg;
455  typedef typename XprType::Rhs RhsArg;
456  typedef evaluator<LhsArg> LhsEvaluator;
457  typedef typename evaluator<RhsArg>::InnerIterator RhsIterator;
458  typedef typename XprType::StorageIndex StorageIndex;
459  typedef typename traits<XprType>::Scalar Scalar;
460 
461  public:
462  class InnerIterator {
463  enum { IsRowMajor = (int(RhsArg::Flags) & RowMajorBit) == RowMajorBit };
464 
465  public:
466  EIGEN_STRONG_INLINE InnerIterator(const sparse_conjunction_evaluator& aEval, Index outer)
467  : m_lhsEval(aEval.m_lhsImpl), m_rhsIter(aEval.m_rhsImpl, outer), m_functor(aEval.m_functor), m_outer(outer) {}
468 
469  EIGEN_STRONG_INLINE InnerIterator& operator++() {
470  ++m_rhsIter;
471  return *this;
472  }
473 
474  EIGEN_STRONG_INLINE Scalar value() const {
475  return m_functor(
476  m_lhsEval.coeff(IsRowMajor ? m_outer : m_rhsIter.index(), IsRowMajor ? m_rhsIter.index() : m_outer),
477  m_rhsIter.value());
478  }
479 
480  EIGEN_STRONG_INLINE StorageIndex index() const { return m_rhsIter.index(); }
481  EIGEN_STRONG_INLINE Index outer() const { return m_rhsIter.outer(); }
482  EIGEN_STRONG_INLINE Index row() const { return m_rhsIter.row(); }
483  EIGEN_STRONG_INLINE Index col() const { return m_rhsIter.col(); }
484 
485  EIGEN_STRONG_INLINE operator bool() const { return m_rhsIter; }
486 
487  protected:
488  const LhsEvaluator& m_lhsEval;
489  RhsIterator m_rhsIter;
490  const BinaryOp& m_functor;
491  const Index m_outer;
492  };
493 
494  enum {
495  CoeffReadCost = int(evaluator<LhsArg>::CoeffReadCost) + int(evaluator<RhsArg>::CoeffReadCost) +
496  int(functor_traits<BinaryOp>::Cost),
497  Flags = XprType::Flags
498  };
499 
500  explicit sparse_conjunction_evaluator(const XprType& xpr)
501  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()) {
502  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
503  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
504  }
505 
506  inline Index nonZerosEstimate() const { return m_rhsImpl.nonZerosEstimate(); }
507 
508  protected:
509  const BinaryOp m_functor;
510  evaluator<LhsArg> m_lhsImpl;
511  evaluator<RhsArg> m_rhsImpl;
512 };
513 
514 // "sparse ^ dense"
515 template <typename XprType>
516 struct sparse_conjunction_evaluator<XprType, IteratorBased, IndexBased> : evaluator_base<XprType> {
517  protected:
518  typedef typename XprType::Functor BinaryOp;
519  typedef typename XprType::Lhs LhsArg;
520  typedef typename XprType::Rhs RhsArg;
521  typedef typename evaluator<LhsArg>::InnerIterator LhsIterator;
522  typedef evaluator<RhsArg> RhsEvaluator;
523  typedef typename XprType::StorageIndex StorageIndex;
524  typedef typename traits<XprType>::Scalar Scalar;
525 
526  public:
527  class InnerIterator {
528  enum { IsRowMajor = (int(LhsArg::Flags) & RowMajorBit) == RowMajorBit };
529 
530  public:
531  EIGEN_STRONG_INLINE InnerIterator(const sparse_conjunction_evaluator& aEval, Index outer)
532  : m_lhsIter(aEval.m_lhsImpl, outer), m_rhsEval(aEval.m_rhsImpl), m_functor(aEval.m_functor), m_outer(outer) {}
533 
534  EIGEN_STRONG_INLINE InnerIterator& operator++() {
535  ++m_lhsIter;
536  return *this;
537  }
538 
539  EIGEN_STRONG_INLINE Scalar value() const {
540  return m_functor(m_lhsIter.value(), m_rhsEval.coeff(IsRowMajor ? m_outer : m_lhsIter.index(),
541  IsRowMajor ? m_lhsIter.index() : m_outer));
542  }
543 
544  EIGEN_STRONG_INLINE StorageIndex index() const { return m_lhsIter.index(); }
545  EIGEN_STRONG_INLINE Index outer() const { return m_lhsIter.outer(); }
546  EIGEN_STRONG_INLINE Index row() const { return m_lhsIter.row(); }
547  EIGEN_STRONG_INLINE Index col() const { return m_lhsIter.col(); }
548 
549  EIGEN_STRONG_INLINE operator bool() const { return m_lhsIter; }
550 
551  protected:
552  LhsIterator m_lhsIter;
553  const evaluator<RhsArg>& m_rhsEval;
554  const BinaryOp& m_functor;
555  const Index m_outer;
556  };
557 
558  enum {
559  CoeffReadCost = int(evaluator<LhsArg>::CoeffReadCost) + int(evaluator<RhsArg>::CoeffReadCost) +
560  int(functor_traits<BinaryOp>::Cost),
561  Flags = XprType::Flags
562  };
563 
564  explicit sparse_conjunction_evaluator(const XprType& xpr)
565  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()) {
566  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
567  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
568  }
569 
570  inline Index nonZerosEstimate() const { return m_lhsImpl.nonZerosEstimate(); }
571 
572  protected:
573  const BinaryOp m_functor;
574  evaluator<LhsArg> m_lhsImpl;
575  evaluator<RhsArg> m_rhsImpl;
576 };
577 
578 template <typename T, typename LhsKind = typename evaluator_traits<typename T::Lhs>::Kind,
579  typename RhsKind = typename evaluator_traits<typename T::Rhs>::Kind,
580  typename LhsScalar = typename traits<typename T::Lhs>::Scalar,
581  typename RhsScalar = typename traits<typename T::Rhs>::Scalar>
582 struct sparse_disjunction_evaluator;
583 
584 // The disjunction "v" evaluator performs a logical "or" or set "union" operation on two input arrays.
585 // Given input arrays 'lhs' and 'rhs' and binary functor 'func',
586 // the sparse destination array 'dst' is evaluated as follows:
587 // if lhs(i,j) and rhs(i,j) are present, dst(i,j) = func(lhs(i,j), rhs(i,j))
588 // if lhs(i,j) is present and rhs(i,j) is null, dst(i,j) = lhs(i,j)
589 // if lhs(i,j) is null and rhs(i,j) is present, dst(i,j) = rhs(i,j)
590 
591 // "sparse v sparse"
592 template <typename XprType>
593 struct sparse_disjunction_evaluator<XprType, IteratorBased, IteratorBased> : evaluator_base<XprType> {
594  protected:
595  typedef typename XprType::Functor BinaryOp;
596  typedef typename XprType::Lhs LhsArg;
597  typedef typename XprType::Rhs RhsArg;
598  typedef typename evaluator<LhsArg>::InnerIterator LhsIterator;
599  typedef typename evaluator<RhsArg>::InnerIterator RhsIterator;
600  typedef typename XprType::StorageIndex StorageIndex;
601  typedef typename traits<XprType>::Scalar Scalar;
602 
603  public:
604  class InnerIterator {
605  public:
606  EIGEN_STRONG_INLINE InnerIterator(const sparse_disjunction_evaluator& aEval, Index outer)
607  : m_lhsIter(aEval.m_lhsImpl, outer),
608  m_rhsIter(aEval.m_rhsImpl, outer),
609  m_functor(aEval.m_functor),
610  m_value(Scalar(0)) {
611  this->operator++();
612  }
613 
614  EIGEN_STRONG_INLINE InnerIterator& operator++() {
615  if (m_lhsIter && m_rhsIter && (m_lhsIter.index() == m_rhsIter.index())) {
616  m_id = m_lhsIter.index();
617  m_value = m_functor(m_lhsIter.value(), m_rhsIter.value());
618  ++m_lhsIter;
619  ++m_rhsIter;
620  } else if (m_lhsIter && (!m_rhsIter || (m_lhsIter.index() < m_rhsIter.index()))) {
621  m_id = m_lhsIter.index();
622  m_value = m_lhsIter.value();
623  ++m_lhsIter;
624  } else if (m_rhsIter && (!m_lhsIter || (m_lhsIter.index() > m_rhsIter.index()))) {
625  m_id = m_rhsIter.index();
626  m_value = m_rhsIter.value();
627  ++m_rhsIter;
628  } else {
629  m_id = -1;
630  }
631  return *this;
632  }
633 
634  EIGEN_STRONG_INLINE Scalar value() const { return m_value; }
635 
636  EIGEN_STRONG_INLINE StorageIndex index() const { return m_id; }
637  EIGEN_STRONG_INLINE Index outer() const { return m_lhsIter.outer(); }
638  EIGEN_STRONG_INLINE Index row() const { return LhsArg::IsRowMajor ? m_lhsIter.row() : index(); }
639  EIGEN_STRONG_INLINE Index col() const { return LhsArg::IsRowMajor ? index() : m_lhsIter.col(); }
640 
641  EIGEN_STRONG_INLINE operator bool() const { return m_id >= 0; }
642 
643  protected:
644  LhsIterator m_lhsIter;
645  RhsIterator m_rhsIter;
646  const BinaryOp& m_functor;
647  Scalar m_value;
648  StorageIndex m_id;
649  };
650 
651  enum {
652  CoeffReadCost = int(evaluator<LhsArg>::CoeffReadCost) + int(evaluator<RhsArg>::CoeffReadCost) +
653  int(functor_traits<BinaryOp>::Cost),
654  Flags = XprType::Flags
655  };
656 
657  explicit sparse_disjunction_evaluator(const XprType& xpr)
658  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()) {
659  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
660  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
661  }
662 
663  inline Index nonZerosEstimate() const { return m_lhsImpl.nonZerosEstimate() + m_rhsImpl.nonZerosEstimate(); }
664 
665  protected:
666  const BinaryOp m_functor;
667  evaluator<LhsArg> m_lhsImpl;
668  evaluator<RhsArg> m_rhsImpl;
669 };
670 
671 // "dense v sparse"
672 template <typename XprType>
673 struct sparse_disjunction_evaluator<XprType, IndexBased, IteratorBased> : evaluator_base<XprType> {
674  protected:
675  typedef typename XprType::Functor BinaryOp;
676  typedef typename XprType::Lhs LhsArg;
677  typedef typename XprType::Rhs RhsArg;
678  typedef evaluator<LhsArg> LhsEvaluator;
679  typedef typename evaluator<RhsArg>::InnerIterator RhsIterator;
680  typedef typename XprType::StorageIndex StorageIndex;
681  typedef typename traits<XprType>::Scalar Scalar;
682 
683  public:
684  class InnerIterator {
685  enum { IsRowMajor = (int(RhsArg::Flags) & RowMajorBit) == RowMajorBit };
686 
687  public:
688  EIGEN_STRONG_INLINE InnerIterator(const sparse_disjunction_evaluator& aEval, Index outer)
689  : m_lhsEval(aEval.m_lhsImpl),
690  m_rhsIter(aEval.m_rhsImpl, outer),
691  m_functor(aEval.m_functor),
692  m_value(0),
693  m_id(-1),
694  m_innerSize(aEval.m_expr.rhs().innerSize()) {
695  this->operator++();
696  }
697 
698  EIGEN_STRONG_INLINE InnerIterator& operator++() {
699  ++m_id;
700  if (m_id < m_innerSize) {
701  Scalar lhsVal = m_lhsEval.coeff(IsRowMajor ? m_rhsIter.outer() : m_id, IsRowMajor ? m_id : m_rhsIter.outer());
702  if (m_rhsIter && m_rhsIter.index() == m_id) {
703  m_value = m_functor(lhsVal, m_rhsIter.value());
704  ++m_rhsIter;
705  } else
706  m_value = lhsVal;
707  }
708 
709  return *this;
710  }
711 
712  EIGEN_STRONG_INLINE Scalar value() const {
713  eigen_internal_assert(m_id < m_innerSize);
714  return m_value;
715  }
716 
717  EIGEN_STRONG_INLINE StorageIndex index() const { return m_id; }
718  EIGEN_STRONG_INLINE Index outer() const { return m_rhsIter.outer(); }
719  EIGEN_STRONG_INLINE Index row() const { return IsRowMajor ? m_rhsIter.outer() : m_id; }
720  EIGEN_STRONG_INLINE Index col() const { return IsRowMajor ? m_id : m_rhsIter.outer(); }
721 
722  EIGEN_STRONG_INLINE operator bool() const { return m_id < m_innerSize; }
723 
724  protected:
725  const evaluator<LhsArg>& m_lhsEval;
726  RhsIterator m_rhsIter;
727  const BinaryOp& m_functor;
728  Scalar m_value;
729  StorageIndex m_id;
730  StorageIndex m_innerSize;
731  };
732 
733  enum {
734  CoeffReadCost = int(evaluator<LhsArg>::CoeffReadCost) + int(evaluator<RhsArg>::CoeffReadCost) +
735  int(functor_traits<BinaryOp>::Cost),
736  Flags = XprType::Flags
737  };
738 
739  explicit sparse_disjunction_evaluator(const XprType& xpr)
740  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()), m_expr(xpr) {
741  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
742  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
743  }
744 
745  inline Index nonZerosEstimate() const { return m_expr.size(); }
746 
747  protected:
748  const BinaryOp m_functor;
749  evaluator<LhsArg> m_lhsImpl;
750  evaluator<RhsArg> m_rhsImpl;
751  const XprType& m_expr;
752 };
753 
754 // "sparse v dense"
755 template <typename XprType>
756 struct sparse_disjunction_evaluator<XprType, IteratorBased, IndexBased> : evaluator_base<XprType> {
757  protected:
758  typedef typename XprType::Functor BinaryOp;
759  typedef typename XprType::Lhs LhsArg;
760  typedef typename XprType::Rhs RhsArg;
761  typedef typename evaluator<LhsArg>::InnerIterator LhsIterator;
762  typedef evaluator<RhsArg> RhsEvaluator;
763  typedef typename XprType::StorageIndex StorageIndex;
764  typedef typename traits<XprType>::Scalar Scalar;
765 
766  public:
767  class InnerIterator {
768  enum { IsRowMajor = (int(LhsArg::Flags) & RowMajorBit) == RowMajorBit };
769 
770  public:
771  EIGEN_STRONG_INLINE InnerIterator(const sparse_disjunction_evaluator& aEval, Index outer)
772  : m_lhsIter(aEval.m_lhsImpl, outer),
773  m_rhsEval(aEval.m_rhsImpl),
774  m_functor(aEval.m_functor),
775  m_value(0),
776  m_id(-1),
777  m_innerSize(aEval.m_expr.lhs().innerSize()) {
778  this->operator++();
779  }
780 
781  EIGEN_STRONG_INLINE InnerIterator& operator++() {
782  ++m_id;
783  if (m_id < m_innerSize) {
784  Scalar rhsVal = m_rhsEval.coeff(IsRowMajor ? m_lhsIter.outer() : m_id, IsRowMajor ? m_id : m_lhsIter.outer());
785  if (m_lhsIter && m_lhsIter.index() == m_id) {
786  m_value = m_functor(m_lhsIter.value(), rhsVal);
787  ++m_lhsIter;
788  } else
789  m_value = rhsVal;
790  }
791 
792  return *this;
793  }
794 
795  EIGEN_STRONG_INLINE Scalar value() const {
796  eigen_internal_assert(m_id < m_innerSize);
797  return m_value;
798  }
799 
800  EIGEN_STRONG_INLINE StorageIndex index() const { return m_id; }
801  EIGEN_STRONG_INLINE Index outer() const { return m_lhsIter.outer(); }
802  EIGEN_STRONG_INLINE Index row() const { return IsRowMajor ? m_lhsIter.outer() : m_id; }
803  EIGEN_STRONG_INLINE Index col() const { return IsRowMajor ? m_id : m_lhsIter.outer(); }
804 
805  EIGEN_STRONG_INLINE operator bool() const { return m_id < m_innerSize; }
806 
807  protected:
808  LhsIterator m_lhsIter;
809  const evaluator<RhsArg>& m_rhsEval;
810  const BinaryOp& m_functor;
811  Scalar m_value;
812  StorageIndex m_id;
813  StorageIndex m_innerSize;
814  };
815 
816  enum {
817  CoeffReadCost = int(evaluator<LhsArg>::CoeffReadCost) + int(evaluator<RhsArg>::CoeffReadCost) +
818  int(functor_traits<BinaryOp>::Cost),
819  Flags = XprType::Flags
820  };
821 
822  explicit sparse_disjunction_evaluator(const XprType& xpr)
823  : m_functor(xpr.functor()), m_lhsImpl(xpr.lhs()), m_rhsImpl(xpr.rhs()), m_expr(xpr) {
824  EIGEN_INTERNAL_CHECK_COST_VALUE(functor_traits<BinaryOp>::Cost);
825  EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
826  }
827 
828  inline Index nonZerosEstimate() const { return m_expr.size(); }
829 
830  protected:
831  const BinaryOp m_functor;
832  evaluator<LhsArg> m_lhsImpl;
833  evaluator<RhsArg> m_rhsImpl;
834  const XprType& m_expr;
835 };
836 
837 // when DupFunc is wrapped with scalar_dup_op, use disjunction evaluator
838 template <typename T1, typename T2, typename DupFunc, typename Lhs, typename Rhs>
839 struct binary_evaluator<CwiseBinaryOp<scalar_disjunction_op<DupFunc, T1, T2>, Lhs, Rhs>, IteratorBased, IteratorBased>
840  : sparse_disjunction_evaluator<CwiseBinaryOp<scalar_disjunction_op<DupFunc, T1, T2>, Lhs, Rhs> > {
841  typedef CwiseBinaryOp<scalar_disjunction_op<DupFunc, T1, T2>, Lhs, Rhs> XprType;
842  typedef sparse_disjunction_evaluator<XprType> Base;
843  explicit binary_evaluator(const XprType& xpr) : Base(xpr) {}
844 };
845 } // namespace internal
846 
847 /***************************************************************************
848  * Implementation of SparseMatrixBase and SparseCwise functions/operators
849  ***************************************************************************/
850 
851 template <typename Derived>
852 template <typename OtherDerived>
853 Derived& SparseMatrixBase<Derived>::operator+=(const EigenBase<OtherDerived>& other) {
854  call_assignment(derived(), other.derived(), internal::add_assign_op<Scalar, typename OtherDerived::Scalar>());
855  return derived();
856 }
857 
858 template <typename Derived>
859 template <typename OtherDerived>
860 Derived& SparseMatrixBase<Derived>::operator-=(const EigenBase<OtherDerived>& other) {
861  call_assignment(derived(), other.derived(), internal::assign_op<Scalar, typename OtherDerived::Scalar>());
862  return derived();
863 }
864 
865 template <typename Derived>
866 template <typename OtherDerived>
867 EIGEN_STRONG_INLINE Derived& SparseMatrixBase<Derived>::operator-=(const SparseMatrixBase<OtherDerived>& other) {
868  return derived() = derived() - other.derived();
869 }
870 
871 template <typename Derived>
872 template <typename OtherDerived>
873 EIGEN_STRONG_INLINE Derived& SparseMatrixBase<Derived>::operator+=(const SparseMatrixBase<OtherDerived>& other) {
874  return derived() = derived() + other.derived();
875 }
876 
877 template <typename Derived>
878 template <typename OtherDerived>
879 Derived& SparseMatrixBase<Derived>::operator+=(const DiagonalBase<OtherDerived>& other) {
880  call_assignment_no_alias(derived(), other.derived(),
881  internal::add_assign_op<Scalar, typename OtherDerived::Scalar>());
882  return derived();
883 }
884 
885 template <typename Derived>
886 template <typename OtherDerived>
887 Derived& SparseMatrixBase<Derived>::operator-=(const DiagonalBase<OtherDerived>& other) {
888  call_assignment_no_alias(derived(), other.derived(),
889  internal::sub_assign_op<Scalar, typename OtherDerived::Scalar>());
890  return derived();
891 }
892 
893 template <typename Derived>
894 template <typename OtherDerived>
895 EIGEN_STRONG_INLINE const typename SparseMatrixBase<Derived>::template CwiseProductDenseReturnType<OtherDerived>::Type
896 SparseMatrixBase<Derived>::cwiseProduct(const MatrixBase<OtherDerived>& other) const {
897  return typename CwiseProductDenseReturnType<OtherDerived>::Type(derived(), other.derived());
898 }
899 
900 template <typename DenseDerived, typename SparseDerived>
901 EIGEN_STRONG_INLINE const
902  CwiseBinaryOp<internal::scalar_sum_op<typename DenseDerived::Scalar, typename SparseDerived::Scalar>,
903  const DenseDerived, const SparseDerived>
904  operator+(const MatrixBase<DenseDerived>& a, const SparseMatrixBase<SparseDerived>& b) {
905  return CwiseBinaryOp<internal::scalar_sum_op<typename DenseDerived::Scalar, typename SparseDerived::Scalar>,
906  const DenseDerived, const SparseDerived>(a.derived(), b.derived());
907 }
908 
909 template <typename SparseDerived, typename DenseDerived>
910 EIGEN_STRONG_INLINE const
911  CwiseBinaryOp<internal::scalar_sum_op<typename SparseDerived::Scalar, typename DenseDerived::Scalar>,
912  const SparseDerived, const DenseDerived>
913  operator+(const SparseMatrixBase<SparseDerived>& a, const MatrixBase<DenseDerived>& b) {
914  return CwiseBinaryOp<internal::scalar_sum_op<typename SparseDerived::Scalar, typename DenseDerived::Scalar>,
915  const SparseDerived, const DenseDerived>(a.derived(), b.derived());
916 }
917 
918 template <typename DenseDerived, typename SparseDerived>
919 EIGEN_STRONG_INLINE const
920  CwiseBinaryOp<internal::scalar_difference_op<typename DenseDerived::Scalar, typename SparseDerived::Scalar>,
921  const DenseDerived, const SparseDerived>
922  operator-(const MatrixBase<DenseDerived>& a, const SparseMatrixBase<SparseDerived>& b) {
923  return CwiseBinaryOp<internal::scalar_difference_op<typename DenseDerived::Scalar, typename SparseDerived::Scalar>,
924  const DenseDerived, const SparseDerived>(a.derived(), b.derived());
925 }
926 
927 template <typename SparseDerived, typename DenseDerived>
928 EIGEN_STRONG_INLINE const
929  CwiseBinaryOp<internal::scalar_difference_op<typename SparseDerived::Scalar, typename DenseDerived::Scalar>,
930  const SparseDerived, const DenseDerived>
931  operator-(const SparseMatrixBase<SparseDerived>& a, const MatrixBase<DenseDerived>& b) {
932  return CwiseBinaryOp<internal::scalar_difference_op<typename SparseDerived::Scalar, typename DenseDerived::Scalar>,
933  const SparseDerived, const DenseDerived>(a.derived(), b.derived());
934 }
935 
936 } // end namespace Eigen
937 
938 #endif // EIGEN_SPARSE_CWISE_BINARY_OP_H
const CwiseBinaryOp< internal::scalar_product_op< Derived ::Scalar, OtherDerived ::Scalar >, const Derived, const OtherDerived > cwiseProduct(const Eigen::SparseMatrixBase< OtherDerived > &other) const
Definition: SparseMatrixBase.h:23
Namespace containing all symbols from the Eigen library.
Definition: B01_Experimental.dox:1
const unsigned int RowMajorBit
Definition: Constants.h:70
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:82