Merge of itasc branch. Project files, scons and cmake should be working. Makefile...
[blender.git] / extern / Eigen2 / Eigen / src / Sparse / SparseCwiseBinaryOp.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra. Eigen itself is part of the KDE project.
3 //
4 // Copyright (C) 2008 Gael Guennebaud <g.gael@free.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_CWISE_BINARY_OP_H
26 #define EIGEN_SPARSE_CWISE_BINARY_OP_H
27
28 // Here we have to handle 3 cases:
29 //  1 - sparse op dense
30 //  2 - dense op sparse
31 //  3 - sparse op sparse
32 // We also need to implement a 4th iterator for:
33 //  4 - dense op dense
34 // Finally, we also need to distinguish between the product and other operations :
35 //                configuration      returned mode
36 //  1 - sparse op dense    product      sparse
37 //                         generic      dense
38 //  2 - dense op sparse    product      sparse
39 //                         generic      dense
40 //  3 - sparse op sparse   product      sparse
41 //                         generic      sparse
42 //  4 - dense op dense     product      dense
43 //                         generic      dense
44
45 template<typename BinaryOp, typename Lhs, typename Rhs>
46 struct ei_traits<SparseCwiseBinaryOp<BinaryOp, Lhs, Rhs> >
47 {
48   typedef typename ei_result_of<
49                      BinaryOp(
50                        typename Lhs::Scalar,
51                        typename Rhs::Scalar
52                      )
53                    >::type Scalar;
54   typedef typename Lhs::Nested LhsNested;
55   typedef typename Rhs::Nested RhsNested;
56   typedef typename ei_unref<LhsNested>::type _LhsNested;
57   typedef typename ei_unref<RhsNested>::type _RhsNested;
58   enum {
59     LhsCoeffReadCost = _LhsNested::CoeffReadCost,
60     RhsCoeffReadCost = _RhsNested::CoeffReadCost,
61     LhsFlags = _LhsNested::Flags,
62     RhsFlags = _RhsNested::Flags,
63     RowsAtCompileTime = Lhs::RowsAtCompileTime,
64     ColsAtCompileTime = Lhs::ColsAtCompileTime,
65     MaxRowsAtCompileTime = Lhs::MaxRowsAtCompileTime,
66     MaxColsAtCompileTime = Lhs::MaxColsAtCompileTime,
67     Flags = (int(LhsFlags) | int(RhsFlags)) & HereditaryBits,
68     CoeffReadCost = LhsCoeffReadCost + RhsCoeffReadCost + ei_functor_traits<BinaryOp>::Cost
69   };
70 };
71
72 template<typename BinaryOp, typename Lhs, typename Rhs>
73 class SparseCwiseBinaryOp : ei_no_assignment_operator,
74   public SparseMatrixBase<SparseCwiseBinaryOp<BinaryOp, Lhs, Rhs> >
75 {
76   public:
77
78     class InnerIterator;
79
80     EIGEN_SPARSE_GENERIC_PUBLIC_INTERFACE(SparseCwiseBinaryOp)
81     typedef typename ei_traits<SparseCwiseBinaryOp>::LhsNested LhsNested;
82     typedef typename ei_traits<SparseCwiseBinaryOp>::RhsNested RhsNested;
83     typedef typename ei_unref<LhsNested>::type _LhsNested;
84     typedef typename ei_unref<RhsNested>::type _RhsNested;
85
86     EIGEN_STRONG_INLINE SparseCwiseBinaryOp(const Lhs& lhs, const Rhs& rhs, const BinaryOp& func = BinaryOp())
87       : m_lhs(lhs), m_rhs(rhs), m_functor(func)
88     {
89       EIGEN_STATIC_ASSERT((_LhsNested::Flags&RowMajorBit)==(_RhsNested::Flags&RowMajorBit),
90         BOTH_MATRICES_MUST_HAVE_THE_SAME_STORAGE_ORDER)
91       EIGEN_STATIC_ASSERT((ei_functor_allows_mixing_real_and_complex<BinaryOp>::ret
92                            ? int(ei_is_same_type<typename Lhs::RealScalar, typename Rhs::RealScalar>::ret)
93                            : int(ei_is_same_type<typename Lhs::Scalar, typename Rhs::Scalar>::ret)),
94         YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
95       // require the sizes to match
96       EIGEN_STATIC_ASSERT_SAME_MATRIX_SIZE(Lhs, Rhs)
97       ei_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
98     }
99
100     EIGEN_STRONG_INLINE int rows() const { return m_lhs.rows(); }
101     EIGEN_STRONG_INLINE int cols() const { return m_lhs.cols(); }
102
103     EIGEN_STRONG_INLINE const _LhsNested& lhs() const { return m_lhs; }
104     EIGEN_STRONG_INLINE const _RhsNested& rhs() const { return m_rhs; }
105     EIGEN_STRONG_INLINE const BinaryOp& functor() const { return m_functor; }
106
107   protected:
108     const LhsNested m_lhs;
109     const RhsNested m_rhs;
110     const BinaryOp m_functor;
111 };
112
113 template<typename BinaryOp, typename Lhs, typename Rhs, typename Derived,
114   int _LhsStorageMode = int(Lhs::Flags) & SparseBit,
115   int _RhsStorageMode = int(Rhs::Flags) & SparseBit>
116 class ei_sparse_cwise_binary_op_inner_iterator_selector;
117
118 template<typename BinaryOp, typename Lhs, typename Rhs>
119 class SparseCwiseBinaryOp<BinaryOp,Lhs,Rhs>::InnerIterator
120   : public ei_sparse_cwise_binary_op_inner_iterator_selector<BinaryOp,Lhs,Rhs, typename SparseCwiseBinaryOp<BinaryOp,Lhs,Rhs>::InnerIterator>
121 {
122   public:
123     typedef ei_sparse_cwise_binary_op_inner_iterator_selector<
124       BinaryOp,Lhs,Rhs, InnerIterator> Base;
125
126     EIGEN_STRONG_INLINE InnerIterator(const SparseCwiseBinaryOp& binOp, int outer)
127       : Base(binOp,outer)
128     {}
129 };
130
131 /***************************************************************************
132 * Implementation of inner-iterators
133 ***************************************************************************/
134
135 // template<typename T> struct ei_func_is_conjunction { enum { ret = false }; };
136 // template<typename T> struct ei_func_is_conjunction<ei_scalar_product_op<T> > { enum { ret = true }; };
137
138 // TODO generalize the ei_scalar_product_op specialization to all conjunctions if any !
139
140 // sparse - sparse  (generic)
141 template<typename BinaryOp, typename Lhs, typename Rhs, typename Derived>
142 class ei_sparse_cwise_binary_op_inner_iterator_selector<BinaryOp, Lhs, Rhs, Derived, IsSparse, IsSparse>
143 {
144     typedef SparseCwiseBinaryOp<BinaryOp, Lhs, Rhs> CwiseBinaryXpr;
145     typedef typename ei_traits<CwiseBinaryXpr>::Scalar Scalar;
146     typedef typename ei_traits<CwiseBinaryXpr>::_LhsNested _LhsNested;
147     typedef typename ei_traits<CwiseBinaryXpr>::_RhsNested _RhsNested;
148     typedef typename _LhsNested::InnerIterator LhsIterator;
149     typedef typename _RhsNested::InnerIterator RhsIterator;
150   public:
151
152     EIGEN_STRONG_INLINE ei_sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, int outer)
153       : m_lhsIter(xpr.lhs(),outer), m_rhsIter(xpr.rhs(),outer), m_functor(xpr.functor())
154     {
155       this->operator++();
156     }
157
158     EIGEN_STRONG_INLINE Derived& operator++()
159     {
160       if (m_lhsIter && m_rhsIter && (m_lhsIter.index() == m_rhsIter.index()))
161       {
162         m_id = m_lhsIter.index();
163         m_value = m_functor(m_lhsIter.value(), m_rhsIter.value());
164         ++m_lhsIter;
165         ++m_rhsIter;
166       }
167       else if (m_lhsIter && (!m_rhsIter || (m_lhsIter.index() < m_rhsIter.index())))
168       {
169         m_id = m_lhsIter.index();
170         m_value = m_functor(m_lhsIter.value(), Scalar(0));
171         ++m_lhsIter;
172       }
173       else if (m_rhsIter && (!m_lhsIter || (m_lhsIter.index() > m_rhsIter.index())))
174       {
175         m_id = m_rhsIter.index();
176         m_value = m_functor(Scalar(0), m_rhsIter.value());
177         ++m_rhsIter;
178       }
179       else
180       {
181         m_id = -1;
182       }
183       return *static_cast<Derived*>(this);
184     }
185
186     EIGEN_STRONG_INLINE Scalar value() const { return m_value; }
187
188     EIGEN_STRONG_INLINE int index() const { return m_id; }
189     EIGEN_STRONG_INLINE int row() const { return m_lhsIter.row(); }
190     EIGEN_STRONG_INLINE int col() const { return m_lhsIter.col(); }
191
192     EIGEN_STRONG_INLINE operator bool() const { return m_id>=0; }
193
194   protected:
195     LhsIterator m_lhsIter;
196     RhsIterator m_rhsIter;
197     const BinaryOp& m_functor;
198     Scalar m_value;
199     int m_id;
200 };
201
202 // sparse - sparse  (product)
203 template<typename T, typename Lhs, typename Rhs, typename Derived>
204 class ei_sparse_cwise_binary_op_inner_iterator_selector<ei_scalar_product_op<T>, Lhs, Rhs, Derived, IsSparse, IsSparse>
205 {
206     typedef ei_scalar_product_op<T> BinaryFunc;
207     typedef SparseCwiseBinaryOp<BinaryFunc, Lhs, Rhs> CwiseBinaryXpr;
208     typedef typename CwiseBinaryXpr::Scalar Scalar;
209     typedef typename ei_traits<CwiseBinaryXpr>::_LhsNested _LhsNested;
210     typedef typename _LhsNested::InnerIterator LhsIterator;
211     typedef typename ei_traits<CwiseBinaryXpr>::_RhsNested _RhsNested;
212     typedef typename _RhsNested::InnerIterator RhsIterator;
213   public:
214
215     EIGEN_STRONG_INLINE ei_sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, int outer)
216       : m_lhsIter(xpr.lhs(),outer), m_rhsIter(xpr.rhs(),outer), m_functor(xpr.functor())
217     {
218       while (m_lhsIter && m_rhsIter && (m_lhsIter.index() != m_rhsIter.index()))
219       {
220         if (m_lhsIter.index() < m_rhsIter.index())
221           ++m_lhsIter;
222         else
223           ++m_rhsIter;
224       }
225     }
226
227     EIGEN_STRONG_INLINE Derived& operator++()
228     {
229       ++m_lhsIter;
230       ++m_rhsIter;
231       while (m_lhsIter && m_rhsIter && (m_lhsIter.index() != m_rhsIter.index()))
232       {
233         if (m_lhsIter.index() < m_rhsIter.index())
234           ++m_lhsIter;
235         else
236           ++m_rhsIter;
237       }
238       return *static_cast<Derived*>(this);
239     }
240
241     EIGEN_STRONG_INLINE Scalar value() const { return m_functor(m_lhsIter.value(), m_rhsIter.value()); }
242
243     EIGEN_STRONG_INLINE int index() const { return m_lhsIter.index(); }
244     EIGEN_STRONG_INLINE int row() const { return m_lhsIter.row(); }
245     EIGEN_STRONG_INLINE int col() const { return m_lhsIter.col(); }
246
247     EIGEN_STRONG_INLINE operator bool() const { return (m_lhsIter && m_rhsIter); }
248
249   protected:
250     LhsIterator m_lhsIter;
251     RhsIterator m_rhsIter;
252     const BinaryFunc& m_functor;
253 };
254
255 // sparse - dense  (product)
256 template<typename T, typename Lhs, typename Rhs, typename Derived>
257 class ei_sparse_cwise_binary_op_inner_iterator_selector<ei_scalar_product_op<T>, Lhs, Rhs, Derived, IsSparse, IsDense>
258 {
259     typedef ei_scalar_product_op<T> BinaryFunc;
260     typedef SparseCwiseBinaryOp<BinaryFunc, Lhs, Rhs> CwiseBinaryXpr;
261     typedef typename CwiseBinaryXpr::Scalar Scalar;
262     typedef typename ei_traits<CwiseBinaryXpr>::_LhsNested _LhsNested;
263     typedef typename ei_traits<CwiseBinaryXpr>::RhsNested RhsNested;
264     typedef typename _LhsNested::InnerIterator LhsIterator;
265     enum { IsRowMajor = (int(Lhs::Flags)&RowMajorBit)==RowMajorBit };
266   public:
267
268     EIGEN_STRONG_INLINE ei_sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, int outer)
269       : m_rhs(xpr.rhs()), m_lhsIter(xpr.lhs(),outer), m_functor(xpr.functor()), m_outer(outer)
270     {}
271
272     EIGEN_STRONG_INLINE Derived& operator++()
273     {
274       ++m_lhsIter;
275       return *static_cast<Derived*>(this);
276     }
277
278     EIGEN_STRONG_INLINE Scalar value() const
279     { return m_functor(m_lhsIter.value(),
280                        m_rhs.coeff(IsRowMajor?m_outer:m_lhsIter.index(),IsRowMajor?m_lhsIter.index():m_outer)); }
281
282     EIGEN_STRONG_INLINE int index() const { return m_lhsIter.index(); }
283     EIGEN_STRONG_INLINE int row() const { return m_lhsIter.row(); }
284     EIGEN_STRONG_INLINE int col() const { return m_lhsIter.col(); }
285
286     EIGEN_STRONG_INLINE operator bool() const { return m_lhsIter; }
287
288   protected:
289     const RhsNested m_rhs;
290     LhsIterator m_lhsIter;
291     const BinaryFunc m_functor;
292     const int m_outer;
293 };
294
295 // sparse - dense  (product)
296 template<typename T, typename Lhs, typename Rhs, typename Derived>
297 class ei_sparse_cwise_binary_op_inner_iterator_selector<ei_scalar_product_op<T>, Lhs, Rhs, Derived, IsDense, IsSparse>
298 {
299     typedef ei_scalar_product_op<T> BinaryFunc;
300     typedef SparseCwiseBinaryOp<BinaryFunc, Lhs, Rhs> CwiseBinaryXpr;
301     typedef typename CwiseBinaryXpr::Scalar Scalar;
302     typedef typename ei_traits<CwiseBinaryXpr>::_RhsNested _RhsNested;
303     typedef typename _RhsNested::InnerIterator RhsIterator;
304     enum { IsRowMajor = (int(Rhs::Flags)&RowMajorBit)==RowMajorBit };
305   public:
306
307     EIGEN_STRONG_INLINE ei_sparse_cwise_binary_op_inner_iterator_selector(const CwiseBinaryXpr& xpr, int outer)
308       : m_xpr(xpr), m_rhsIter(xpr.rhs(),outer), m_functor(xpr.functor()), m_outer(outer)
309     {}
310
311     EIGEN_STRONG_INLINE Derived& operator++()
312     {
313       ++m_rhsIter;
314       return *static_cast<Derived*>(this);
315     }
316
317     EIGEN_STRONG_INLINE Scalar value() const
318     { return m_functor(m_xpr.lhs().coeff(IsRowMajor?m_outer:m_rhsIter.index(),IsRowMajor?m_rhsIter.index():m_outer), m_rhsIter.value()); }
319
320     EIGEN_STRONG_INLINE int index() const { return m_rhsIter.index(); }
321     EIGEN_STRONG_INLINE int row() const { return m_rhsIter.row(); }
322     EIGEN_STRONG_INLINE int col() const { return m_rhsIter.col(); }
323
324     EIGEN_STRONG_INLINE operator bool() const { return m_rhsIter; }
325
326   protected:
327     const CwiseBinaryXpr& m_xpr;
328     RhsIterator m_rhsIter;
329     const BinaryFunc& m_functor;
330     const int m_outer;
331 };
332
333
334 /***************************************************************************
335 * Implementation of SparseMatrixBase and SparseCwise functions/operators
336 ***************************************************************************/
337
338 template<typename Derived>
339 template<typename OtherDerived>
340 EIGEN_STRONG_INLINE const SparseCwiseBinaryOp<ei_scalar_difference_op<typename ei_traits<Derived>::Scalar>,
341                                  Derived, OtherDerived>
342 SparseMatrixBase<Derived>::operator-(const SparseMatrixBase<OtherDerived> &other) const
343 {
344   return SparseCwiseBinaryOp<ei_scalar_difference_op<Scalar>,
345                        Derived, OtherDerived>(derived(), other.derived());
346 }
347
348 template<typename Derived>
349 template<typename OtherDerived>
350 EIGEN_STRONG_INLINE Derived &
351 SparseMatrixBase<Derived>::operator-=(const SparseMatrixBase<OtherDerived> &other)
352 {
353   return *this = derived() - other.derived();
354 }
355
356 template<typename Derived>
357 template<typename OtherDerived>
358 EIGEN_STRONG_INLINE const SparseCwiseBinaryOp<ei_scalar_sum_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
359 SparseMatrixBase<Derived>::operator+(const SparseMatrixBase<OtherDerived> &other) const
360 {
361   return SparseCwiseBinaryOp<ei_scalar_sum_op<Scalar>, Derived, OtherDerived>(derived(), other.derived());
362 }
363
364 template<typename Derived>
365 template<typename OtherDerived>
366 EIGEN_STRONG_INLINE Derived &
367 SparseMatrixBase<Derived>::operator+=(const SparseMatrixBase<OtherDerived>& other)
368 {
369   return *this = derived() + other.derived();
370 }
371
372 template<typename ExpressionType>
373 template<typename OtherDerived>
374 EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE
375 SparseCwise<ExpressionType>::operator*(const SparseMatrixBase<OtherDerived> &other) const
376 {
377   return EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE(_expression(), other.derived());
378 }
379
380 template<typename ExpressionType>
381 template<typename OtherDerived>
382 EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE
383 SparseCwise<ExpressionType>::operator*(const MatrixBase<OtherDerived> &other) const
384 {
385   return EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE(_expression(), other.derived());
386 }
387
388 // template<typename ExpressionType>
389 // template<typename OtherDerived>
390 // EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_quotient_op)
391 // SparseCwise<ExpressionType>::operator/(const SparseMatrixBase<OtherDerived> &other) const
392 // {
393 //   return EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_quotient_op)(_expression(), other.derived());
394 // }
395 //
396 // template<typename ExpressionType>
397 // template<typename OtherDerived>
398 // EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_quotient_op)
399 // SparseCwise<ExpressionType>::operator/(const MatrixBase<OtherDerived> &other) const
400 // {
401 //   return EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_quotient_op)(_expression(), other.derived());
402 // }
403
404 template<typename ExpressionType>
405 template<typename OtherDerived>
406 inline ExpressionType& SparseCwise<ExpressionType>::operator*=(const SparseMatrixBase<OtherDerived> &other)
407 {
408   return m_matrix.const_cast_derived() = _expression() * other.derived();
409 }
410
411 // template<typename ExpressionType>
412 // template<typename OtherDerived>
413 // inline ExpressionType& SparseCwise<ExpressionType>::operator/=(const SparseMatrixBase<OtherDerived> &other)
414 // {
415 //   return m_matrix.const_cast_derived() = *this / other;
416 // }
417
418 template<typename ExpressionType>
419 template<typename OtherDerived>
420 EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_min_op)
421 SparseCwise<ExpressionType>::min(const SparseMatrixBase<OtherDerived> &other) const
422 {
423   return EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_min_op)(_expression(), other.derived());
424 }
425
426 template<typename ExpressionType>
427 template<typename OtherDerived>
428 EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_max_op)
429 SparseCwise<ExpressionType>::max(const SparseMatrixBase<OtherDerived> &other) const
430 {
431   return EIGEN_SPARSE_CWISE_BINOP_RETURN_TYPE(ei_scalar_max_op)(_expression(), other.derived());
432 }
433
434 // template<typename Derived>
435 // template<typename CustomBinaryOp, typename OtherDerived>
436 // EIGEN_STRONG_INLINE const CwiseBinaryOp<CustomBinaryOp, Derived, OtherDerived>
437 // SparseMatrixBase<Derived>::binaryExpr(const SparseMatrixBase<OtherDerived> &other, const CustomBinaryOp& func) const
438 // {
439 //   return CwiseBinaryOp<CustomBinaryOp, Derived, OtherDerived>(derived(), other.derived(), func);
440 // }
441
442 #endif // EIGEN_SPARSE_CWISE_BINARY_OP_H