SparseMatrixBase.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-2011 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_SPARSEMATRIXBASE_H
26 #define EIGEN_SPARSEMATRIXBASE_H
27 
28 namespace Eigen {
29 
41 template<typename Derived> class SparseMatrixBase : public EigenBase<Derived>
42 {
43  public:
44 
45  typedef typename internal::traits<Derived>::Scalar Scalar;
46  typedef typename internal::packet_traits<Scalar>::type PacketScalar;
47  typedef typename internal::traits<Derived>::StorageKind StorageKind;
48  typedef typename internal::traits<Derived>::Index Index;
49  typedef typename internal::add_const_on_value_type_if_arithmetic<
50  typename internal::packet_traits<Scalar>::type
52 
55 
56  template<typename OtherDerived>
57  Derived& operator=(const EigenBase<OtherDerived> &other)
58  {
59  other.derived().evalTo(derived());
60  return derived();
61  }
62 
63  enum {
64 
65  RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime,
71  ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime,
78  SizeAtCompileTime = (internal::size_at_compile_time<internal::traits<Derived>::RowsAtCompileTime,
79  internal::traits<Derived>::ColsAtCompileTime>::ret),
86 
87  MaxSizeAtCompileTime = (internal::size_at_compile_time<MaxRowsAtCompileTime,
88  MaxColsAtCompileTime>::ret),
89 
96  Flags = internal::traits<Derived>::Flags,
101  CoeffReadCost = internal::traits<Derived>::CoeffReadCost,
107 
108  #ifndef EIGEN_PARSED_BY_DOXYGEN
109  _HasDirectAccess = (int(Flags)&DirectAccessBit) ? 1 : 0 // workaround sunCC
110  #endif
111  };
112 
118 
119 
121 
122 
123 #ifndef EIGEN_PARSED_BY_DOXYGEN
130  typedef typename NumTraits<Scalar>::Real RealScalar;
131 
134  typedef typename internal::conditional<_HasDirectAccess, const Scalar&, Scalar>::type CoeffReturnType;
135 
138 
142 
143  inline const Derived& derived() const { return *static_cast<const Derived*>(this); }
144  inline Derived& derived() { return *static_cast<Derived*>(this); }
145  inline Derived& const_cast_derived() const
146  { return *static_cast<Derived*>(const_cast<SparseMatrixBase*>(this)); }
147 #endif // not EIGEN_PARSED_BY_DOXYGEN
148 
149 #define EIGEN_CURRENT_STORAGE_BASE_CLASS Eigen::SparseMatrixBase
150 # include "../plugins/CommonCwiseUnaryOps.h"
151 # include "../plugins/CommonCwiseBinaryOps.h"
152 # include "../plugins/MatrixCwiseUnaryOps.h"
153 # include "../plugins/MatrixCwiseBinaryOps.h"
154 # ifdef EIGEN_SPARSEMATRIXBASE_PLUGIN
155 # include EIGEN_SPARSEMATRIXBASE_PLUGIN
156 # endif
157 # undef EIGEN_CURRENT_STORAGE_BASE_CLASS
158 #undef EIGEN_CURRENT_STORAGE_BASE_CLASS
159 
160 
162  inline Index rows() const { return derived().rows(); }
164  inline Index cols() const { return derived().cols(); }
167  inline Index size() const { return rows() * cols(); }
170  inline Index nonZeros() const { return derived().nonZeros(); }
175  inline bool isVector() const { return rows()==1 || cols()==1; }
178  Index outerSize() const { return (int(Flags)&RowMajorBit) ? this->rows() : this->cols(); }
181  Index innerSize() const { return (int(Flags)&RowMajorBit) ? this->cols() : this->rows(); }
183  bool isRValue() const { return m_isRValue; }
184  Derived& markAsRValue() { m_isRValue = true; return derived(); }
185 
186  SparseMatrixBase() : m_isRValue(false) { /* TODO check flags */ }
187 
188 
189  template<typename OtherDerived>
190  Derived& operator=(const ReturnByValue<OtherDerived>& other)
191  {
192  other.evalTo(derived());
193  return derived();
194  }
195 
196 
197  template<typename OtherDerived>
198  inline Derived& operator=(const SparseMatrixBase<OtherDerived>& other)
199  {
200  return assign(other.derived());
201  }
202 
203  inline Derived& operator=(const Derived& other)
204  {
205 // if (other.isRValue())
206 // derived().swap(other.const_cast_derived());
207 // else
208  return assign(other.derived());
209  }
210 
211  protected:
212 
213  template<typename OtherDerived>
214  inline Derived& assign(const OtherDerived& other)
215  {
216  const bool transpose = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit);
217  const Index outerSize = (int(OtherDerived::Flags) & RowMajorBit) ? other.rows() : other.cols();
218  if ((!transpose) && other.isRValue())
219  {
220  // eval without temporary
221  derived().resize(other.rows(), other.cols());
222  derived().setZero();
223  derived().reserve((std::max)(this->rows(),this->cols())*2);
224  for (Index j=0; j<outerSize; ++j)
225  {
226  derived().startVec(j);
227  for (typename OtherDerived::InnerIterator it(other, j); it; ++it)
228  {
229  Scalar v = it.value();
230  derived().insertBackByOuterInner(j,it.index()) = v;
231  }
232  }
233  derived().finalize();
234  }
235  else
236  {
237  assignGeneric(other);
238  }
239  return derived();
240  }
241 
242  template<typename OtherDerived>
243  inline void assignGeneric(const OtherDerived& other)
244  {
245  //const bool transpose = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit);
246  eigen_assert(( ((internal::traits<Derived>::SupportedAccessPatterns&OuterRandomAccessPattern)==OuterRandomAccessPattern) ||
247  (!((Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit)))) &&
248  "the transpose operation is supposed to be handled in SparseMatrix::operator=");
249 
250  enum { Flip = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit) };
251 
252  const Index outerSize = other.outerSize();
253  //typedef typename internal::conditional<transpose, LinkedVectorMatrix<Scalar,Flags&RowMajorBit>, Derived>::type TempType;
254  // thanks to shallow copies, we always eval to a tempary
255  Derived temp(other.rows(), other.cols());
256 
257  temp.reserve((std::max)(this->rows(),this->cols())*2);
258  for (Index j=0; j<outerSize; ++j)
259  {
260  temp.startVec(j);
261  for (typename OtherDerived::InnerIterator it(other.derived(), j); it; ++it)
262  {
263  Scalar v = it.value();
264  temp.insertBackByOuterInner(Flip?it.index():j,Flip?j:it.index()) = v;
265  }
266  }
267  temp.finalize();
268 
269  derived() = temp.markAsRValue();
270  }
271 
272  public:
273 
274  template<typename Lhs, typename Rhs>
275  inline Derived& operator=(const SparseSparseProduct<Lhs,Rhs>& product);
276 
277  friend std::ostream & operator << (std::ostream & s, const SparseMatrixBase& m)
278  {
279  typedef typename Derived::Nested Nested;
280  typedef typename internal::remove_all<Nested>::type NestedCleaned;
281 
282  if (Flags&RowMajorBit)
283  {
284  const Nested nm(m.derived());
285  for (Index row=0; row<nm.outerSize(); ++row)
286  {
287  Index col = 0;
288  for (typename NestedCleaned::InnerIterator it(nm.derived(), row); it; ++it)
289  {
290  for ( ; col<it.index(); ++col)
291  s << "0 ";
292  s << it.value() << " ";
293  ++col;
294  }
295  for ( ; col<m.cols(); ++col)
296  s << "0 ";
297  s << std::endl;
298  }
299  }
300  else
301  {
302  const Nested nm(m.derived());
303  if (m.cols() == 1) {
304  Index row = 0;
305  for (typename NestedCleaned::InnerIterator it(nm.derived(), 0); it; ++it)
306  {
307  for ( ; row<it.index(); ++row)
308  s << "0" << std::endl;
309  s << it.value() << std::endl;
310  ++row;
311  }
312  for ( ; row<m.rows(); ++row)
313  s << "0" << std::endl;
314  }
315  else
316  {
318  s << static_cast<const SparseMatrixBase<SparseMatrix<Scalar, RowMajorBit> >&>(trans);
319  }
320  }
321  return s;
322  }
323 
324  template<typename OtherDerived>
325  Derived& operator+=(const SparseMatrixBase<OtherDerived>& other);
326  template<typename OtherDerived>
327  Derived& operator-=(const SparseMatrixBase<OtherDerived>& other);
328 
329  Derived& operator*=(const Scalar& other);
330  Derived& operator/=(const Scalar& other);
331 
332  #define EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE \
333  CwiseBinaryOp< \
334  internal::scalar_product_op< \
335  typename internal::scalar_product_traits< \
336  typename internal::traits<Derived>::Scalar, \
337  typename internal::traits<OtherDerived>::Scalar \
338  >::ReturnType \
339  >, \
340  Derived, \
341  OtherDerived \
342  >
343 
344  template<typename OtherDerived>
346  cwiseProduct(const MatrixBase<OtherDerived> &other) const;
347 
348  // sparse * sparse
349  template<typename OtherDerived>
351  operator*(const SparseMatrixBase<OtherDerived> &other) const;
352 
353  // sparse * diagonal
354  template<typename OtherDerived>
355  const SparseDiagonalProduct<Derived,OtherDerived>
356  operator*(const DiagonalBase<OtherDerived> &other) const;
357 
358  // diagonal * sparse
359  template<typename OtherDerived> friend
360  const SparseDiagonalProduct<OtherDerived,Derived>
362  { return SparseDiagonalProduct<OtherDerived,Derived>(lhs.derived(), rhs.derived()); }
363 
365  template<typename OtherDerived> friend
367  operator*(const MatrixBase<OtherDerived>& lhs, const Derived& rhs)
368  { return typename DenseSparseProductReturnType<OtherDerived,Derived>::Type(lhs.derived(),rhs); }
369 
371  template<typename OtherDerived>
373  operator*(const MatrixBase<OtherDerived> &other) const;
374 
377  {
379  }
380 
381  template<typename OtherDerived>
382  Derived& operator*=(const SparseMatrixBase<OtherDerived>& other);
383 
384  #ifdef EIGEN2_SUPPORT
385  // deprecated
386  template<typename OtherDerived>
388  solveTriangular(const MatrixBase<OtherDerived>& other) const;
389 
390  // deprecated
391  template<typename OtherDerived>
392  void solveTriangularInPlace(MatrixBase<OtherDerived>& other) const;
393  #endif // EIGEN2_SUPPORT
394 
395  template<int Mode>
397 
398  template<unsigned int UpLo> inline const SparseSelfAdjointView<Derived, UpLo> selfadjointView() const;
399  template<unsigned int UpLo> inline SparseSelfAdjointView<Derived, UpLo> selfadjointView();
400 
401  template<typename OtherDerived> Scalar dot(const MatrixBase<OtherDerived>& other) const;
402  template<typename OtherDerived> Scalar dot(const SparseMatrixBase<OtherDerived>& other) const;
403  RealScalar squaredNorm() const;
404  RealScalar norm() const;
405 
407  const Transpose<const Derived> transpose() const { return derived(); }
408  const AdjointReturnType adjoint() const { return transpose(); }
409 
410  // sub-vector
417 
418  // set of sub-vectors
423 
430 
432  template<typename DenseDerived>
434  {
435  dst.setZero();
436  for (Index j=0; j<outerSize(); ++j)
437  for (typename Derived::InnerIterator i(derived(),j); i; ++i)
438  dst.coeffRef(i.row(),i.col()) = i.value();
439  }
440 
442  {
443  return derived();
444  }
445 
446  template<typename OtherDerived>
448  RealScalar prec = NumTraits<Scalar>::dummy_precision()) const
449  { return toDense().isApprox(other.toDense(),prec); }
450 
451  template<typename OtherDerived>
453  RealScalar prec = NumTraits<Scalar>::dummy_precision()) const
454  { return toDense().isApprox(other,prec); }
455 
461  inline const typename internal::eval<Derived>::type eval() const
462  { return typename internal::eval<Derived>::type(derived()); }
463 
464  Scalar sum() const;
465 
466  protected:
467 
469 };
470 
471 } // end namespace Eigen
472 
473 #endif // EIGEN_SPARSEMATRIXBASE_H