blob: 7220beea4f02e277e786f9087af496827ff29c2b [file] [log] [blame]
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +03001// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
Gael Guennebaud3c634462014-07-01 13:27:35 +02004// Copyright (C) 2011-2014 Gael Guennebaud <gael.guennebaud@inria.fr>
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +03005// Copyright (C) 2010 Daniel Lowengrub <lowdanie@gmail.com>
6//
Benoit Jacob69124cf2012-07-13 14:42:47 -04007// This Source Code Form is subject to the terms of the Mozilla
8// Public License v. 2.0. If a copy of the MPL was not distributed
9// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +030010
11#ifndef EIGEN_SPARSEVIEW_H
12#define EIGEN_SPARSEVIEW_H
13
Antonio Sánchez6e4d5d42023-08-21 16:25:22 +000014// IWYU pragma: private
Rasmus Munk Larsend7d0bf82021-09-10 19:12:26 +000015#include "./InternalHeaderCheck.h"
16
Tobias Woodf38e16c2023-11-29 11:12:48 +000017namespace Eigen {
Jitse Niesen3c412182012-04-15 11:06:28 +010018
Benoit Jacob47160402010-10-25 10:15:22 -040019namespace internal {
20
Tobias Woodf38e16c2023-11-29 11:12:48 +000021template <typename MatrixType>
22struct traits<SparseView<MatrixType> > : traits<MatrixType> {
Christoph Hertzberge8cdbed2014-12-04 22:48:53 +010023 typedef typename MatrixType::StorageIndex StorageIndex;
Gael Guennebaud99e4afd2010-06-24 22:48:48 +020024 typedef Sparse StorageKind;
Tobias Woodf38e16c2023-11-29 11:12:48 +000025 enum { Flags = int(traits<MatrixType>::Flags) & (RowMajorBit) };
Gael Guennebaud99e4afd2010-06-24 22:48:48 +020026};
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +030027
Tobias Woodf38e16c2023-11-29 11:12:48 +000028} // end namespace internal
Benoit Jacob47160402010-10-25 10:15:22 -040029
Gael Guennebaud831fffe2017-01-06 18:01:29 +010030/** \ingroup SparseCore_Module
Tobias Woodf38e16c2023-11-29 11:12:48 +000031 * \class SparseView
32 *
33 * \brief Expression of a dense or sparse matrix with zero or too small values removed
34 *
35 * \tparam MatrixType the type of the object of which we are removing the small entries
36 *
37 * This class represents an expression of a given dense or sparse matrix with
38 * entries smaller than \c reference * \c epsilon are removed.
39 * It is the return type of MatrixBase::sparseView() and SparseMatrixBase::pruned()
40 * and most of the time this is the only way it is used.
41 *
42 * \sa MatrixBase::sparseView(), SparseMatrixBase::pruned()
43 */
44template <typename MatrixType>
45class SparseView : public SparseMatrixBase<SparseView<MatrixType> > {
Gael Guennebaud99e4afd2010-06-24 22:48:48 +020046 typedef typename MatrixType::Nested MatrixTypeNested;
Erik Schultheis421cbf02022-03-16 16:43:40 +000047 typedef internal::remove_all_t<MatrixTypeNested> MatrixTypeNested_;
Tobias Woodf38e16c2023-11-29 11:12:48 +000048 typedef SparseMatrixBase<SparseView> Base;
49
50 public:
Benoit Jacob2d65f5d2010-06-14 09:06:27 -040051 EIGEN_SPARSE_PUBLIC_INTERFACE(SparseView)
Erik Schultheis421cbf02022-03-16 16:43:40 +000052 typedef internal::remove_all_t<MatrixType> NestedExpression;
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +030053
Gael Guennebaud9a2447b2015-06-09 09:11:12 +020054 explicit SparseView(const MatrixType& mat, const Scalar& reference = Scalar(0),
Tobias Woodf38e16c2023-11-29 11:12:48 +000055 const RealScalar& epsilon = NumTraits<Scalar>::dummy_precision())
56 : m_matrix(mat), m_reference(reference), m_epsilon(epsilon) {}
Gael Guennebaudfa6d36e2010-07-22 15:57:01 +020057
Gael Guennebaudfc202ba2015-02-13 18:57:41 +010058 inline Index rows() const { return m_matrix.rows(); }
59 inline Index cols() const { return m_matrix.cols(); }
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +030060
Gael Guennebaudfc202ba2015-02-13 18:57:41 +010061 inline Index innerSize() const { return m_matrix.innerSize(); }
62 inline Index outerSize() const { return m_matrix.outerSize(); }
Tobias Woodf38e16c2023-11-29 11:12:48 +000063
Gael Guennebaudf0648f82014-06-26 13:52:19 +020064 /** \returns the nested expression */
Tobias Woodf38e16c2023-11-29 11:12:48 +000065 const internal::remove_all_t<MatrixTypeNested>& nestedExpression() const { return m_matrix; }
66
Gael Guennebaudf0648f82014-06-26 13:52:19 +020067 Scalar reference() const { return m_reference; }
68 RealScalar epsilon() const { return m_epsilon; }
Tobias Woodf38e16c2023-11-29 11:12:48 +000069
70 protected:
Gael Guennebaudfe85b7e2012-02-03 23:18:26 +010071 MatrixTypeNested m_matrix;
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +030072 Scalar m_reference;
Gael Guennebaudf0648f82014-06-26 13:52:19 +020073 RealScalar m_epsilon;
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +030074};
75
Gael Guennebaudf0648f82014-06-26 13:52:19 +020076namespace internal {
77
78// TODO find a way to unify the two following variants
79// This is tricky because implementing an inner iterator on top of an IndexBased evaluator is
80// not easy because the evaluators do not expose the sizes of the underlying expression.
Gael Guennebaudf0648f82014-06-26 13:52:19 +020081
Tobias Woodf38e16c2023-11-29 11:12:48 +000082template <typename ArgType>
83struct unary_evaluator<SparseView<ArgType>, IteratorBased> : public evaluator_base<SparseView<ArgType> > {
84 typedef typename evaluator<ArgType>::InnerIterator EvalIterator;
Gael Guennebaudf0648f82014-06-26 13:52:19 +020085
Tobias Woodf38e16c2023-11-29 11:12:48 +000086 public:
87 typedef SparseView<ArgType> XprType;
Gael Guennebaudf0648f82014-06-26 13:52:19 +020088
Tobias Woodf38e16c2023-11-29 11:12:48 +000089 class InnerIterator : public EvalIterator {
90 protected:
Gael Guennebaudf0648f82014-06-26 13:52:19 +020091 typedef typename XprType::Scalar Scalar;
Gael Guennebaudf0648f82014-06-26 13:52:19 +020092
Tobias Woodf38e16c2023-11-29 11:12:48 +000093 public:
94 EIGEN_STRONG_INLINE InnerIterator(const unary_evaluator& sve, Index outer)
95 : EvalIterator(sve.m_argImpl, outer), m_view(sve.m_view) {
96 incrementToNonZero();
97 }
Gael Guennebaudf0648f82014-06-26 13:52:19 +020098
Tobias Woodf38e16c2023-11-29 11:12:48 +000099 EIGEN_STRONG_INLINE InnerIterator& operator++() {
100 EvalIterator::operator++();
101 incrementToNonZero();
102 return *this;
103 }
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200104
Tobias Woodf38e16c2023-11-29 11:12:48 +0000105 using EvalIterator::value;
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200106
Tobias Woodf38e16c2023-11-29 11:12:48 +0000107 protected:
108 const XprType& m_view;
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200109
Tobias Woodf38e16c2023-11-29 11:12:48 +0000110 private:
111 void incrementToNonZero() {
112 while ((bool(*this)) && internal::isMuchSmallerThan(value(), m_view.reference(), m_view.epsilon())) {
113 EvalIterator::operator++();
114 }
115 }
116 };
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200117
Tobias Woodf38e16c2023-11-29 11:12:48 +0000118 enum { CoeffReadCost = evaluator<ArgType>::CoeffReadCost, Flags = XprType::Flags };
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200119
Tobias Woodf38e16c2023-11-29 11:12:48 +0000120 explicit unary_evaluator(const XprType& xpr) : m_argImpl(xpr.nestedExpression()), m_view(xpr) {}
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200121
Tobias Woodf38e16c2023-11-29 11:12:48 +0000122 protected:
123 evaluator<ArgType> m_argImpl;
124 const XprType& m_view;
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200125};
126
Tobias Woodf38e16c2023-11-29 11:12:48 +0000127template <typename ArgType>
128struct unary_evaluator<SparseView<ArgType>, IndexBased> : public evaluator_base<SparseView<ArgType> > {
129 public:
130 typedef SparseView<ArgType> XprType;
131
132 protected:
133 enum { IsRowMajor = (XprType::Flags & RowMajorBit) == RowMajorBit };
134 typedef typename XprType::Scalar Scalar;
135 typedef typename XprType::StorageIndex StorageIndex;
136
137 public:
138 class InnerIterator {
139 public:
140 EIGEN_STRONG_INLINE InnerIterator(const unary_evaluator& sve, Index outer)
141 : m_sve(sve), m_inner(0), m_outer(outer), m_end(sve.m_view.innerSize()) {
142 incrementToNonZero();
143 }
144
145 EIGEN_STRONG_INLINE InnerIterator& operator++() {
146 m_inner++;
147 incrementToNonZero();
148 return *this;
149 }
150
151 EIGEN_STRONG_INLINE Scalar value() const {
152 return (IsRowMajor) ? m_sve.m_argImpl.coeff(m_outer, m_inner) : m_sve.m_argImpl.coeff(m_inner, m_outer);
153 }
154
155 EIGEN_STRONG_INLINE StorageIndex index() const { return m_inner; }
156 inline Index row() const { return IsRowMajor ? m_outer : index(); }
157 inline Index col() const { return IsRowMajor ? index() : m_outer; }
158
159 EIGEN_STRONG_INLINE operator bool() const { return m_inner < m_end && m_inner >= 0; }
160
161 protected:
162 const unary_evaluator& m_sve;
163 Index m_inner;
164 const Index m_outer;
165 const Index m_end;
166
167 private:
168 void incrementToNonZero() {
169 while ((bool(*this)) && internal::isMuchSmallerThan(value(), m_sve.m_view.reference(), m_sve.m_view.epsilon())) {
170 m_inner++;
171 }
172 }
173 };
174
175 enum { CoeffReadCost = evaluator<ArgType>::CoeffReadCost, Flags = XprType::Flags };
176
177 explicit unary_evaluator(const XprType& xpr) : m_argImpl(xpr.nestedExpression()), m_view(xpr) {}
178
179 protected:
180 evaluator<ArgType> m_argImpl;
181 const XprType& m_view;
182};
183
184} // end namespace internal
Gael Guennebaudf0648f82014-06-26 13:52:19 +0200185
Gael Guennebaud831fffe2017-01-06 18:01:29 +0100186/** \ingroup SparseCore_Module
Tobias Woodf38e16c2023-11-29 11:12:48 +0000187 *
188 * \returns a sparse expression of the dense expression \c *this with values smaller than
189 * \a reference * \a epsilon removed.
190 *
191 * This method is typically used when prototyping to convert a quickly assembled dense Matrix \c D to a SparseMatrix \c
192 * S: \code MatrixXd D(n,m); SparseMatrix<double> S; S = D.sparseView(); // suppress numerical zeros (exact)
193 * S = D.sparseView(reference);
194 * S = D.sparseView(reference,epsilon);
195 * \endcode
196 * where \a reference is a meaningful non zero reference value,
197 * and \a epsilon is a tolerance factor defaulting to NumTraits<Scalar>::dummy_precision().
198 *
199 * \sa SparseMatrixBase::pruned(), class SparseView */
200template <typename Derived>
Gael Guennebaud746d2db2014-07-01 13:18:56 +0200201const SparseView<Derived> MatrixBase<Derived>::sparseView(const Scalar& reference,
Tobias Woodf38e16c2023-11-29 11:12:48 +0000202 const typename NumTraits<Scalar>::Real& epsilon) const {
Gael Guennebaud746d2db2014-07-01 13:18:56 +0200203 return SparseView<Derived>(derived(), reference, epsilon);
204}
205
206/** \returns an expression of \c *this with values smaller than
Tobias Woodf38e16c2023-11-29 11:12:48 +0000207 * \a reference * \a epsilon removed.
208 *
209 * This method is typically used in conjunction with the product of two sparse matrices
210 * to automatically prune the smallest values as follows:
211 * \code
212 * C = (A*B).pruned(); // suppress numerical zeros (exact)
213 * C = (A*B).pruned(ref);
214 * C = (A*B).pruned(ref,epsilon);
215 * \endcode
216 * where \c ref is a meaningful non zero reference value.
217 * */
218template <typename Derived>
219const SparseView<Derived> SparseMatrixBase<Derived>::pruned(const Scalar& reference, const RealScalar& epsilon) const {
Gael Guennebaud746d2db2014-07-01 13:18:56 +0200220 return SparseView<Derived>(derived(), reference, epsilon);
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +0300221}
222
Tobias Woodf38e16c2023-11-29 11:12:48 +0000223} // end namespace Eigen
Jitse Niesen3c412182012-04-15 11:06:28 +0100224
Daniel Lowengrubdcd39a92010-06-14 02:16:46 +0300225#endif