DynamicSparseMatrix.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008-2009 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_DYNAMIC_SPARSEMATRIX_H
26 #define EIGEN_DYNAMIC_SPARSEMATRIX_H
27 
28 namespace Eigen {
29 
50 namespace internal {
51 template<typename _Scalar, int _Options, typename _Index>
52 struct traits<DynamicSparseMatrix<_Scalar, _Options, _Index> >
53 {
54  typedef _Scalar Scalar;
55  typedef _Index Index;
56  typedef Sparse StorageKind;
57  typedef MatrixXpr XprKind;
58  enum {
59  RowsAtCompileTime = Dynamic,
60  ColsAtCompileTime = Dynamic,
61  MaxRowsAtCompileTime = Dynamic,
62  MaxColsAtCompileTime = Dynamic,
63  Flags = _Options | NestByRefBit | LvalueBit,
64  CoeffReadCost = NumTraits<Scalar>::ReadCost,
65  SupportedAccessPatterns = OuterRandomAccessPattern
66  };
67 };
68 }
69 
70 template<typename _Scalar, int _Options, typename _Index>
72  : public SparseMatrixBase<DynamicSparseMatrix<_Scalar, _Options, _Index> >
73 {
74  public:
75  EIGEN_SPARSE_PUBLIC_INTERFACE(DynamicSparseMatrix)
76  // FIXME: why are these operator already alvailable ???
77  // EIGEN_SPARSE_INHERIT_ASSIGNMENT_OPERATOR(DynamicSparseMatrix, +=)
78  // EIGEN_SPARSE_INHERIT_ASSIGNMENT_OPERATOR(DynamicSparseMatrix, -=)
79  typedef MappedSparseMatrix<Scalar,Flags> Map;
80  using Base::IsRowMajor;
81  using Base::operator=;
82  enum {
83  Options = _Options
84  };
85 
86  protected:
87 
89 
90  Index m_innerSize;
91  std::vector<internal::CompressedStorage<Scalar,Index> > m_data;
92 
93  public:
94 
95  inline Index rows() const { return IsRowMajor ? outerSize() : m_innerSize; }
96  inline Index cols() const { return IsRowMajor ? m_innerSize : outerSize(); }
97  inline Index innerSize() const { return m_innerSize; }
98  inline Index outerSize() const { return static_cast<Index>(m_data.size()); }
99  inline Index innerNonZeros(Index j) const { return m_data[j].size(); }
100 
101  std::vector<internal::CompressedStorage<Scalar,Index> >& _data() { return m_data; }
102  const std::vector<internal::CompressedStorage<Scalar,Index> >& _data() const { return m_data; }
103 
107  inline Scalar coeff(Index row, Index col) const
108  {
109  const Index outer = IsRowMajor ? row : col;
110  const Index inner = IsRowMajor ? col : row;
111  return m_data[outer].at(inner);
112  }
113 
118  inline Scalar& coeffRef(Index row, Index col)
119  {
120  const Index outer = IsRowMajor ? row : col;
121  const Index inner = IsRowMajor ? col : row;
122  return m_data[outer].atWithInsertion(inner);
123  }
124 
125  class InnerIterator;
126  class ReverseInnerIterator;
127 
128  void setZero()
129  {
130  for (Index j=0; j<outerSize(); ++j)
131  m_data[j].clear();
132  }
133 
135  Index nonZeros() const
136  {
137  Index res = 0;
138  for (Index j=0; j<outerSize(); ++j)
139  res += static_cast<Index>(m_data[j].size());
140  return res;
141  }
142 
143 
144 
145  void reserve(Index reserveSize = 1000)
146  {
147  if (outerSize()>0)
148  {
149  Index reserveSizePerVector = (std::max)(reserveSize/outerSize(),Index(4));
150  for (Index j=0; j<outerSize(); ++j)
151  {
152  m_data[j].reserve(reserveSizePerVector);
153  }
154  }
155  }
156 
158  inline void startVec(Index /*outer*/) {}
159 
165  inline Scalar& insertBack(Index row, Index col)
166  {
167  return insertBackByOuterInner(IsRowMajor?row:col, IsRowMajor?col:row);
168  }
169 
171  inline Scalar& insertBackByOuterInner(Index outer, Index inner)
172  {
173  eigen_assert(outer<Index(m_data.size()) && inner<m_innerSize && "out of range");
174  eigen_assert(((m_data[outer].size()==0) || (m_data[outer].index(m_data[outer].size()-1)<inner))
175  && "wrong sorted insertion");
176  m_data[outer].append(0, inner);
177  return m_data[outer].value(m_data[outer].size()-1);
178  }
179 
180  inline Scalar& insert(Index row, Index col)
181  {
182  const Index outer = IsRowMajor ? row : col;
183  const Index inner = IsRowMajor ? col : row;
184 
185  Index startId = 0;
186  Index id = static_cast<Index>(m_data[outer].size()) - 1;
187  m_data[outer].resize(id+2,1);
188 
189  while ( (id >= startId) && (m_data[outer].index(id) > inner) )
190  {
191  m_data[outer].index(id+1) = m_data[outer].index(id);
192  m_data[outer].value(id+1) = m_data[outer].value(id);
193  --id;
194  }
195  m_data[outer].index(id+1) = inner;
196  m_data[outer].value(id+1) = 0;
197  return m_data[outer].value(id+1);
198  }
199 
201  inline void finalize() {}
202 
204  void prune(Scalar reference, RealScalar epsilon = NumTraits<RealScalar>::dummy_precision())
205  {
206  for (Index j=0; j<outerSize(); ++j)
207  m_data[j].prune(reference,epsilon);
208  }
209 
212  void resize(Index rows, Index cols)
213  {
214  const Index outerSize = IsRowMajor ? rows : cols;
215  m_innerSize = IsRowMajor ? cols : rows;
216  setZero();
217  if (Index(m_data.size()) != outerSize)
218  {
219  m_data.resize(outerSize);
220  }
221  }
222 
223  void resizeAndKeepData(Index rows, Index cols)
224  {
225  const Index outerSize = IsRowMajor ? rows : cols;
226  const Index innerSize = IsRowMajor ? cols : rows;
227  if (m_innerSize>innerSize)
228  {
229  // remove all coefficients with innerCoord>=innerSize
230  // TODO
231  //std::cerr << "not implemented yet\n";
232  exit(2);
233  }
234  if (m_data.size() != outerSize)
235  {
236  m_data.resize(outerSize);
237  }
238  }
239 
241  EIGEN_DEPRECATED inline DynamicSparseMatrix()
242  : m_innerSize(0), m_data(0)
243  {
244  eigen_assert(innerSize()==0 && outerSize()==0);
245  }
246 
248  EIGEN_DEPRECATED inline DynamicSparseMatrix(Index rows, Index cols)
249  : m_innerSize(0)
250  {
251  resize(rows, cols);
252  }
253 
255  template<typename OtherDerived>
256  EIGEN_DEPRECATED explicit inline DynamicSparseMatrix(const SparseMatrixBase<OtherDerived>& other)
257  : m_innerSize(0)
258  {
259  Base::operator=(other.derived());
260  }
261 
262  inline DynamicSparseMatrix(const DynamicSparseMatrix& other)
263  : Base(), m_innerSize(0)
264  {
265  *this = other.derived();
266  }
267 
268  inline void swap(DynamicSparseMatrix& other)
269  {
270  //EIGEN_DBG_SPARSE(std::cout << "SparseMatrix:: swap\n");
271  std::swap(m_innerSize, other.m_innerSize);
272  //std::swap(m_outerSize, other.m_outerSize);
273  m_data.swap(other.m_data);
274  }
275 
276  inline DynamicSparseMatrix& operator=(const DynamicSparseMatrix& other)
277  {
278  if (other.isRValue())
279  {
280  swap(other.const_cast_derived());
281  }
282  else
283  {
284  resize(other.rows(), other.cols());
285  m_data = other.m_data;
286  }
287  return *this;
288  }
289 
292 
293  public:
294 
297  EIGEN_DEPRECATED void startFill(Index reserveSize = 1000)
298  {
299  setZero();
300  reserve(reserveSize);
301  }
302 
312  EIGEN_DEPRECATED Scalar& fill(Index row, Index col)
313  {
314  const Index outer = IsRowMajor ? row : col;
315  const Index inner = IsRowMajor ? col : row;
316  return insertBack(outer,inner);
317  }
318 
324  EIGEN_DEPRECATED Scalar& fillrand(Index row, Index col)
325  {
326  return insert(row,col);
327  }
328 
331  EIGEN_DEPRECATED void endFill() {}
332 
333 # ifdef EIGEN_DYNAMICSPARSEMATRIX_PLUGIN
334 # include EIGEN_DYNAMICSPARSEMATRIX_PLUGIN
335 # endif
336  };
337 
338 template<typename Scalar, int _Options, typename _Index>
339 class DynamicSparseMatrix<Scalar,_Options,_Index>::InnerIterator : public SparseVector<Scalar,_Options,_Index>::InnerIterator
340 {
341  typedef typename SparseVector<Scalar,_Options,_Index>::InnerIterator Base;
342  public:
343  InnerIterator(const DynamicSparseMatrix& mat, Index outer)
344  : Base(mat.m_data[outer]), m_outer(outer)
345  {}
346 
347  inline Index row() const { return IsRowMajor ? m_outer : Base::index(); }
348  inline Index col() const { return IsRowMajor ? Base::index() : m_outer; }
349 
350  protected:
351  const Index m_outer;
352 };
353 
354 template<typename Scalar, int _Options, typename _Index>
355 class DynamicSparseMatrix<Scalar,_Options,_Index>::ReverseInnerIterator : public SparseVector<Scalar,_Options,_Index>::ReverseInnerIterator
356 {
357  typedef typename SparseVector<Scalar,_Options,_Index>::ReverseInnerIterator Base;
358  public:
359  ReverseInnerIterator(const DynamicSparseMatrix& mat, Index outer)
360  : Base(mat.m_data[outer]), m_outer(outer)
361  {}
362 
363  inline Index row() const { return IsRowMajor ? m_outer : Base::index(); }
364  inline Index col() const { return IsRowMajor ? Base::index() : m_outer; }
365 
366  protected:
367  const Index m_outer;
368 };
369 
370 } // end namespace Eigen
371 
372 #endif // EIGEN_DYNAMIC_SPARSEMATRIX_H