aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen
diff options
context:
space:
mode:
authorGravatar Gael Guennebaud <g.gael@free.fr>2008-10-21 13:35:04 +0000
committerGravatar Gael Guennebaud <g.gael@free.fr>2008-10-21 13:35:04 +0000
commitcf0f82ecbe1b92ec44e8fe34f65d6183059c9491 (patch)
treea4d6b680561765f192fc384100c68071e86648ab /Eigen
parent9e02e42ff6757da1c96518a24c913a01250564ee (diff)
sparse module:
- remove some useless stuff => let's focus on a single sparse matrix format - finalize the new RandomSetter
Diffstat (limited to 'Eigen')
-rw-r--r--Eigen/Sparse10
-rw-r--r--Eigen/src/Sparse/HashMatrix.h166
-rw-r--r--Eigen/src/Sparse/LinkedVectorMatrix.h317
-rw-r--r--Eigen/src/Sparse/RandomSetter.h153
-rw-r--r--Eigen/src/Sparse/SparseMatrixBase.h7
-rw-r--r--Eigen/src/Sparse/SparseSetter.h138
6 files changed, 149 insertions, 642 deletions
diff --git a/Eigen/Sparse b/Eigen/Sparse
index 2364c6fe4..54b15b8ae 100644
--- a/Eigen/Sparse
+++ b/Eigen/Sparse
@@ -8,6 +8,10 @@
#include <cstring>
#include <algorithm>
+#ifdef EIGEN_GOOGLEHASH_SUPPORT
+ #include <google/dense_hash_map>
+#endif
+
#ifdef EIGEN_CHOLMOD_SUPPORT
extern "C" {
#include "cholmod.h"
@@ -70,10 +74,10 @@ namespace Eigen {
#include "src/Sparse/RandomSetter.h"
#include "src/Sparse/SparseBlock.h"
#include "src/Sparse/SparseMatrix.h"
-#include "src/Sparse/HashMatrix.h"
-#include "src/Sparse/LinkedVectorMatrix.h"
+//#include "src/Sparse/HashMatrix.h"
+//#include "src/Sparse/LinkedVectorMatrix.h"
#include "src/Sparse/CoreIterators.h"
-#include "src/Sparse/SparseSetter.h"
+//#include "src/Sparse/SparseSetter.h"
#include "src/Sparse/SparseProduct.h"
#include "src/Sparse/TriangularSolver.h"
#include "src/Sparse/SparseLLT.h"
diff --git a/Eigen/src/Sparse/HashMatrix.h b/Eigen/src/Sparse/HashMatrix.h
deleted file mode 100644
index 1d6d33cce..000000000
--- a/Eigen/src/Sparse/HashMatrix.h
+++ /dev/null
@@ -1,166 +0,0 @@
-// This file is part of Eigen, a lightweight C++ template library
-// for linear algebra. Eigen itself is part of the KDE project.
-//
-// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
-//
-// Eigen is free software; you can redistribute it and/or
-// modify it under the terms of the GNU Lesser General Public
-// License as published by the Free Software Foundation; either
-// version 3 of the License, or (at your option) any later version.
-//
-// Alternatively, you can redistribute it and/or
-// modify it under the terms of the GNU General Public License as
-// published by the Free Software Foundation; either version 2 of
-// the License, or (at your option) any later version.
-//
-// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
-// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
-// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
-// GNU General Public License for more details.
-//
-// You should have received a copy of the GNU Lesser General Public
-// License and a copy of the GNU General Public License along with
-// Eigen. If not, see <http://www.gnu.org/licenses/>.
-
-#ifndef EIGEN_HASHMATRIX_H
-#define EIGEN_HASHMATRIX_H
-
-template<typename _Scalar, int _Flags>
-struct ei_traits<HashMatrix<_Scalar, _Flags> >
-{
- typedef _Scalar Scalar;
- enum {
- RowsAtCompileTime = Dynamic,
- ColsAtCompileTime = Dynamic,
- MaxRowsAtCompileTime = Dynamic,
- MaxColsAtCompileTime = Dynamic,
- Flags = SparseBit | _Flags,
- CoeffReadCost = NumTraits<Scalar>::ReadCost,
- SupportedAccessPatterns = RandomAccessPattern
- };
-};
-
-// TODO reimplement this class using custom linked lists
-template<typename _Scalar, int _Flags>
-class HashMatrix
- : public SparseMatrixBase<HashMatrix<_Scalar, _Flags> >
-{
- public:
- EIGEN_GENERIC_PUBLIC_INTERFACE(HashMatrix)
- class InnerIterator;
- protected:
-
- typedef typename std::map<int, Scalar>::iterator MapIterator;
- typedef typename std::map<int, Scalar>::const_iterator ConstMapIterator;
-
- public:
- inline int rows() const { return m_innerSize; }
- inline int cols() const { return m_data.size(); }
-
- inline const Scalar& coeff(int row, int col) const
- {
- const MapIterator it = m_data[col].find(row);
- if (it!=m_data[col].end())
- return Scalar(0);
- return it->second;
- }
-
- inline Scalar& coeffRef(int row, int col)
- {
- return m_data[col][row];
- }
-
- public:
-
- inline void startFill(int /*reserveSize = 1000 --- currently unused, don't generate a warning*/) {}
-
- inline Scalar& fill(int row, int col) { return coeffRef(row, col); }
-
- inline void endFill() {}
-
- ~HashMatrix()
- {}
-
- inline void shallowCopy(const HashMatrix& other)
- {
- EIGEN_DBG_SPARSE(std::cout << "HashMatrix:: shallowCopy\n");
- // FIXME implement a true shallow copy !!
- resize(other.rows(), other.cols());
- for (int j=0; j<this->outerSize(); ++j)
- m_data[j] = other.m_data[j];
- }
-
- void resize(int _rows, int _cols)
- {
- if (cols() != _cols)
- {
- m_data.resize(_cols);
- }
- m_innerSize = _rows;
- }
-
- inline HashMatrix(int rows, int cols)
- : m_innerSize(0)
- {
- resize(rows, cols);
- }
-
- template<typename OtherDerived>
- inline HashMatrix(const MatrixBase<OtherDerived>& other)
- : m_innerSize(0)
- {
- *this = other.derived();
- }
-
- inline HashMatrix& operator=(const HashMatrix& other)
- {
- if (other.isRValue())
- {
- shallowCopy(other);
- }
- else
- {
- resize(other.rows(), other.cols());
- for (int col=0; col<cols(); ++col)
- m_data[col] = other.m_data[col];
- }
- return *this;
- }
-
- template<typename OtherDerived>
- inline HashMatrix& operator=(const MatrixBase<OtherDerived>& other)
- {
- return SparseMatrixBase<HashMatrix>::operator=(other);
- }
-
- protected:
-
- std::vector<std::map<int, Scalar> > m_data;
- int m_innerSize;
-
-};
-
-template<typename Scalar, int _Flags>
-class HashMatrix<Scalar,_Flags>::InnerIterator
-{
- public:
-
- InnerIterator(const HashMatrix& mat, int col)
- : m_matrix(mat), m_it(mat.m_data[col].begin()), m_end(mat.m_data[col].end())
- {}
-
- InnerIterator& operator++() { m_it++; return *this; }
-
- Scalar value() { return m_it->second; }
-
- int index() const { return m_it->first; }
-
- operator bool() const { return m_it!=m_end; }
-
- protected:
- const HashMatrix& m_matrix;
- typename HashMatrix::ConstMapIterator m_it;
- typename HashMatrix::ConstMapIterator m_end;
-};
-
-#endif // EIGEN_HASHMATRIX_H
diff --git a/Eigen/src/Sparse/LinkedVectorMatrix.h b/Eigen/src/Sparse/LinkedVectorMatrix.h
deleted file mode 100644
index 4e5a4862e..000000000
--- a/Eigen/src/Sparse/LinkedVectorMatrix.h
+++ /dev/null
@@ -1,317 +0,0 @@
-// This file is part of Eigen, a lightweight C++ template library
-// for linear algebra. Eigen itself is part of the KDE project.
-//
-// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
-//
-// Eigen is free software; you can redistribute it and/or
-// modify it under the terms of the GNU Lesser General Public
-// License as published by the Free Software Foundation; either
-// version 3 of the License, or (at your option) any later version.
-//
-// Alternatively, you can redistribute it and/or
-// modify it under the terms of the GNU General Public License as
-// published by the Free Software Foundation; either version 2 of
-// the License, or (at your option) any later version.
-//
-// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
-// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
-// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
-// GNU General Public License for more details.
-//
-// You should have received a copy of the GNU Lesser General Public
-// License and a copy of the GNU General Public License along with
-// Eigen. If not, see <http://www.gnu.org/licenses/>.
-
-#ifndef EIGEN_LINKEDVECTORMATRIX_H
-#define EIGEN_LINKEDVECTORMATRIX_H
-
-template<typename _Scalar, int _Flags>
-struct ei_traits<LinkedVectorMatrix<_Scalar,_Flags> >
-{
- typedef _Scalar Scalar;
- enum {
- RowsAtCompileTime = Dynamic,
- ColsAtCompileTime = Dynamic,
- MaxRowsAtCompileTime = Dynamic,
- MaxColsAtCompileTime = Dynamic,
- Flags = SparseBit | _Flags,
- CoeffReadCost = NumTraits<Scalar>::ReadCost,
- SupportedAccessPatterns = InnerCoherentAccessPattern
- };
-};
-
-template<typename Element, int ChunkSize = 8>
-struct LinkedVectorChunk
-{
- LinkedVectorChunk() : next(0), prev(0), size(0) {}
- Element data[ChunkSize];
- LinkedVectorChunk* next;
- LinkedVectorChunk* prev;
- int size;
- bool isFull() const { return size==ChunkSize; }
-};
-
-template<typename _Scalar, int _Flags>
-class LinkedVectorMatrix
- : public SparseMatrixBase<LinkedVectorMatrix<_Scalar,_Flags> >
-{
- public:
- EIGEN_GENERIC_PUBLIC_INTERFACE(LinkedVectorMatrix)
- class InnerIterator;
- protected:
-
- enum {
- RowMajor = Flags&RowMajorBit ? 1 : 0
- };
-
- struct ValueIndex
- {
- ValueIndex() : value(0), index(0) {}
- ValueIndex(Scalar v, int i) : value(v), index(i) {}
- Scalar value;
- int index;
- };
- typedef LinkedVectorChunk<ValueIndex,8> VectorChunk;
-
- inline int find(VectorChunk** _el, int id)
- {
- VectorChunk* el = *_el;
- while (el && el->data[el->size-1].index<id)
- el = el->next;
- *_el = el;
- if (el)
- {
- // binary search
- int maxI = el->size-1;
- int minI = 0;
- int i = el->size/2;
- const ValueIndex* data = el->data;
- while (data[i].index!=id)
- {
- if (data[i].index<id)
- {
- minI = i+1;
- i = (maxI + minI)+2;
- }
- else
- {
- maxI = i-1;
- i = (maxI + minI)+2;
- }
- if (minI>=maxI)
- return -1;
- }
- if (data[i].index==id)
- return i;
- }
- return -1;
- }
-
- public:
- inline int rows() const { return RowMajor ? m_data.size() : m_innerSize; }
- inline int cols() const { return RowMajor ? m_innerSize : m_data.size(); }
-
- inline const Scalar& coeff(int row, int col) const
- {
- const int outer = RowMajor ? row : col;
- const int inner = RowMajor ? col : row;
-
- VectorChunk* el = m_data[outer];
- int id = find(&el, inner);
- if (id<0)
- return Scalar(0);
- return el->data[id].value;
- }
-
- inline Scalar& coeffRef(int row, int col)
- {
- const int outer = RowMajor ? row : col;
- const int inner = RowMajor ? col : row;
-
- VectorChunk* el = m_data[outer];
- int id = find(&el, inner);
- ei_assert(id>=0);
-// if (id<0)
-// return Scalar(0);
- return el->data[id].value;
- }
-
- public:
-
- inline void startFill(int reserveSize = 1000)
- {
- clear();
- for (unsigned int i=0; i<m_data.size(); ++i)
- m_ends[i] = m_data[i] = 0;
- }
-
- inline Scalar& fill(int row, int col)
- {
- const int outer = RowMajor ? row : col;
- const int inner = RowMajor ? col : row;
-// std::cout << " ll fill " << outer << "," << inner << "\n";
- if (m_ends[outer]==0)
- {
- m_data[outer] = m_ends[outer] = new VectorChunk();
- }
- else
- {
- ei_assert(m_ends[outer]->data[m_ends[outer]->size-1].index < inner);
- if (m_ends[outer]->isFull())
- {
-
- VectorChunk* el = new VectorChunk();
- m_ends[outer]->next = el;
- el->prev = m_ends[outer];
- m_ends[outer] = el;
- }
- }
- m_ends[outer]->data[m_ends[outer]->size].index = inner;
- return m_ends[outer]->data[m_ends[outer]->size++].value;
- }
-
- inline void endFill() { }
-
- void printDbg()
- {
- for (int j=0; j<m_data.size(); ++j)
- {
- VectorChunk* el = m_data[j];
- while (el)
- {
- for (int i=0; i<el->size; ++i)
- std::cout << j << "," << el->data[i].index << " = " << el->data[i].value << "\n";
- el = el->next;
- }
- }
- for (int j=0; j<m_data.size(); ++j)
- {
- InnerIterator it(*this,j);
- while (it)
- {
- std::cout << j << "," << it.index() << " = " << it.value() << "\n";
- ++it;
- }
- }
- }
-
- ~LinkedVectorMatrix()
- {
- clear();
- }
-
- void clear()
- {
- for (unsigned int i=0; i<m_data.size(); ++i)
- {
- VectorChunk* el = m_data[i];
- while (el)
- {
- VectorChunk* tmp = el;
- el = el->next;
- delete tmp;
- }
- }
- }
-
- void resize(int rows, int cols)
- {
- const int outers = RowMajor ? rows : cols;
- const int inners = RowMajor ? cols : rows;
-
- if (this->outerSize() != outers)
- {
- clear();
- m_data.resize(outers);
- m_ends.resize(outers);
- for (unsigned int i=0; i<m_data.size(); ++i)
- m_ends[i] = m_data[i] = 0;
- }
- m_innerSize = inners;
- }
-
- inline LinkedVectorMatrix(int rows, int cols)
- : m_innerSize(0)
- {
- resize(rows, cols);
- }
-
- template<typename OtherDerived>
- inline LinkedVectorMatrix(const MatrixBase<OtherDerived>& other)
- : m_innerSize(0)
- {
- *this = other.derived();
- }
-
- inline void swap(LinkedVectorMatrix& other)
- {
- EIGEN_DBG_SPARSE(std::cout << "LinkedVectorMatrix:: swap\n");
- resize(other.rows(), other.cols());
- m_data.swap(other.m_data);
- m_ends.swap(other.m_ends);
- }
-
- inline LinkedVectorMatrix& operator=(const LinkedVectorMatrix& other)
- {
- if (other.isRValue())
- {
- swap(other.const_cast_derived());
- }
- else
- {
- // TODO implement a specialized deep copy here
- return operator=<LinkedVectorMatrix>(other);
- }
- return *this;
- }
-
- template<typename OtherDerived>
- inline LinkedVectorMatrix& operator=(const MatrixBase<OtherDerived>& other)
- {
- return SparseMatrixBase<LinkedVectorMatrix>::operator=(other.derived());
- }
-
- protected:
-
- // outer vector of inner linked vector chunks
- std::vector<VectorChunk*> m_data;
- // stores a reference to the last vector chunk for efficient filling
- std::vector<VectorChunk*> m_ends;
- int m_innerSize;
-
-};
-
-
-template<typename Scalar, int _Flags>
-class LinkedVectorMatrix<Scalar,_Flags>::InnerIterator
-{
- public:
-
- InnerIterator(const LinkedVectorMatrix& mat, int col)
- : m_matrix(mat), m_el(mat.m_data[col]), m_it(0)
- {}
-
- InnerIterator& operator++()
- {
- m_it++;
- if (m_it>=m_el->size)
- {
- m_el = m_el->next;
- m_it = 0;
- }
- return *this;
- }
-
- Scalar value() { return m_el->data[m_it].value; }
-
- int index() const { return m_el->data[m_it].index; }
-
- operator bool() const { return m_el && (m_el->next || m_it<m_el->size); }
-
- protected:
- const LinkedVectorMatrix& m_matrix;
- VectorChunk* m_el;
- int m_it;
-};
-
-#endif // EIGEN_LINKEDVECTORMATRIX_H
diff --git a/Eigen/src/Sparse/RandomSetter.h b/Eigen/src/Sparse/RandomSetter.h
index 8606e03c6..35bc9daee 100644
--- a/Eigen/src/Sparse/RandomSetter.h
+++ b/Eigen/src/Sparse/RandomSetter.h
@@ -29,6 +29,9 @@ template<typename Scalar> struct StdMapTraits
{
typedef int KeyType;
typedef std::map<KeyType,Scalar> Type;
+ enum {
+ IsSorted = 1
+ };
static void setInvalidKey(Type&, const KeyType&) {}
};
@@ -38,6 +41,9 @@ template<typename Scalar> struct GnuHashMapTraits
{
typedef int KeyType;
typedef __gnu_cxx::hash_map<KeyType,Scalar> Type;
+ enum {
+ IsSorted = 0
+ };
static void setInvalidKey(Type&, const KeyType&) {}
};
@@ -48,6 +54,9 @@ template<typename Scalar> struct GoogleDenseHashMapTraits
{
typedef int KeyType;
typedef google::dense_hash_map<KeyType,Scalar> Type;
+ enum {
+ IsSorted = 0
+ };
static void setInvalidKey(Type& map, const KeyType& k)
{ map.set_empty_key(k); }
@@ -59,6 +68,9 @@ template<typename Scalar> struct GoogleSparseHashMapTraits
{
typedef int KeyType;
typedef google::sparse_hash_map<KeyType,Scalar> Type;
+ enum {
+ IsSorted = 0
+ };
static void setInvalidKey(Type&, const KeyType&) {}
};
@@ -66,10 +78,35 @@ template<typename Scalar> struct GoogleSparseHashMapTraits
/** \class RandomSetter
*
+ * Typical usage:
+ * \code
+ * SparseMatrix<double> m(rows,cols);
+ * {
+ * RandomSetter<SparseMatrix<double> > w(m);
+ * // don't use m but w instead with read/write random access to the coefficients:
+ * for(;;)
+ * w(rand(),rand()) = rand;
+ * }
+ * // when w is deleted, the data are copied back to m
+ * // and m is ready to use.
+ * \endcode
+ *
+ * \note for performance and memory consumption reasons it is highly recommended to use
+ * Google's hash library. To do so you have two options:
+ * - include <google/dense_hash_map> yourself \b before Eigen/Sparse header
+ * - define EIGEN_GOOGLEHASH_SUPPORT
+ * In the later case the inclusion of <google/dense_hash_map> is made for you.
*/
template<typename SparseMatrixType,
- template <typename T> class HashMapTraits = StdMapTraits,
- int OuterPacketBits = 6>
+ template <typename T> class MapTraits =
+#if defined _DENSE_HASH_MAP_H_
+ GoogleDenseHashMapTraits
+#elif defined _HASH_MAP
+ GnuHashMapTraits
+#else
+ StdMapTraits
+#endif
+ ,int OuterPacketBits = 6>
class RandomSetter
{
typedef typename ei_traits<SparseMatrixType>::Scalar Scalar;
@@ -78,11 +115,13 @@ class RandomSetter
ScalarWrapper() : value(0) {}
Scalar value;
};
- typedef typename HashMapTraits<ScalarWrapper>::KeyType KeyType;
- typedef typename HashMapTraits<ScalarWrapper>::Type HashMapType;
+ typedef typename MapTraits<ScalarWrapper>::KeyType KeyType;
+ typedef typename MapTraits<ScalarWrapper>::Type HashMapType;
static const int OuterPacketMask = (1 << OuterPacketBits) - 1;
enum {
- RowMajor = SparseMatrixType::Flags & RowMajorBit
+ SwapStorage = 1 - MapTraits<ScalarWrapper>::IsSorted,
+ TargetRowMajor = (SparseMatrixType::Flags & RowMajorBit) ? 1 : 0,
+ SetterRowMajor = SwapStorage ? 1-TargetRowMajor : TargetRowMajor
};
public:
@@ -90,31 +129,114 @@ class RandomSetter
inline RandomSetter(SparseMatrixType& target)
: mp_target(&target)
{
- m_outerPackets = target.outerSize() >> OuterPacketBits;
- if (target.outerSize()&OuterPacketMask)
+ const int outerSize = SwapStorage ? target.innerSize() : target.outerSize();
+ const int innerSize = SwapStorage ? target.outerSize() : target.innerSize();
+ m_outerPackets = outerSize >> OuterPacketBits;
+ if (outerSize&OuterPacketMask)
m_outerPackets += 1;
m_hashmaps = new HashMapType[m_outerPackets];
- KeyType ik = (1<<OuterPacketBits)*mp_target->innerSize()+1;
+ // compute number of bits needed to store inner indices
+ int aux = innerSize - 1;
+ m_keyBitsOffset = 0;
+ while (aux)
+ {
+ m_keyBitsOffset++;
+ aux = aux >> 1;
+ }
+ KeyType ik = (1<<(OuterPacketBits+m_keyBitsOffset));
for (int k=0; k<m_outerPackets; ++k)
- HashMapTraits<ScalarWrapper>::setInvalidKey(m_hashmaps[k],ik);
+ MapTraits<ScalarWrapper>::setInvalidKey(m_hashmaps[k],ik);
+
+ // insert current coeffs
+ for (int j=0; j<mp_target->outerSize(); ++j)
+ for (typename SparseMatrixType::InnerIterator it(*mp_target,j); it; ++it)
+ (*this)(TargetRowMajor?j:it.index(), TargetRowMajor?it.index():j) = it.value();
}
~RandomSetter()
{
+ KeyType keyBitsMask = (1<<m_keyBitsOffset)-1;
+ if (!SwapStorage) // also means the map is sorted
+ {
+ mp_target->startFill(nonZeros());
+ for (int k=0; k<m_outerPackets; ++k)
+ {
+ const int outerOffset = (1<<OuterPacketBits) * k;
+ typename HashMapType::iterator end = m_hashmaps[k].end();
+ for (typename HashMapType::iterator it = m_hashmaps[k].begin(); it!=end; ++it)
+ {
+ const int outer = (it->first >> m_keyBitsOffset) + outerOffset;
+ const int inner = it->first & keyBitsMask;
+ mp_target->fill(TargetRowMajor ? outer : inner, TargetRowMajor ? inner : outer) = it->second.value;
+ }
+ }
+ mp_target->endFill();
+ }
+ else
+ {
+ VectorXi positions(mp_target->outerSize());
+ positions.setZero();
+ // pass 1
+ for (int k=0; k<m_outerPackets; ++k)
+ {
+ typename HashMapType::iterator end = m_hashmaps[k].end();
+ for (typename HashMapType::iterator it = m_hashmaps[k].begin(); it!=end; ++it)
+ {
+ const int outer = it->first & keyBitsMask;
+ positions[outer]++;
+ }
+ }
+ // prefix sum
+ int count = 0;
+ for (int j=0; j<mp_target->outerSize(); ++j)
+ {
+ int tmp = positions[j];
+ mp_target->_outerIndexPtr()[j] = count;
+ positions[j] = count;
+ count += tmp;
+ }
+ mp_target->_outerIndexPtr()[mp_target->outerSize()] = count;
+ mp_target->resizeNonZeros(count);
+ // pass 2
+ for (int k=0; k<m_outerPackets; ++k)
+ {
+ const int outerOffset = (1<<OuterPacketBits) * k;
+ typename HashMapType::iterator end = m_hashmaps[k].end();
+ for (typename HashMapType::iterator it = m_hashmaps[k].begin(); it!=end; ++it)
+ {
+ const int inner = (it->first >> m_keyBitsOffset) + outerOffset;
+ const int outer = it->first & keyBitsMask;
+ // sorted insertion
+ // Note that we have to deal with at most 2^OuterPacketBits unsorted coefficients,
+ // moreover those 2^OuterPacketBits coeffs are likely to be sparse, an so only a
+ // small fraction of them have to be sorted, whence the following simple procedure:
+ int posStart = mp_target->_outerIndexPtr()[outer];
+ int i = (positions[outer]++) - 1;
+ while ( (i >= posStart) && (mp_target->_innerIndexPtr()[i] > inner) )
+ {
+ mp_target->_valuePtr()[i+1] = mp_target->_valuePtr()[i];
+ mp_target->_innerIndexPtr()[i+1] = mp_target->_innerIndexPtr()[i];
+ --i;
+ }
+ mp_target->_innerIndexPtr()[i+1] = inner;
+ mp_target->_valuePtr()[i+1] = it->second.value;
+ }
+ }
+ }
delete[] m_hashmaps;
}
Scalar& operator() (int row, int col)
{
- const int outer = RowMajor ? row : col;
- const int inner = RowMajor ? col : row;
- const int outerMajor = outer >> OuterPacketBits;
- const int outerMinor = outer & OuterPacketMask;
- const KeyType key = inner + outerMinor * mp_target->innerSize();
-
+ const int outer = SetterRowMajor ? row : col;
+ const int inner = SetterRowMajor ? col : row;
+ const int outerMajor = outer >> OuterPacketBits; // index of the packet/map
+ const int outerMinor = outer & OuterPacketMask; // index of the inner vector in the packet
+ const KeyType key = (KeyType(outerMinor)<<m_keyBitsOffset) | inner;
return m_hashmaps[outerMajor][key].value;
}
+ // might be slow
int nonZeros() const
{
int nz = 0;
@@ -129,6 +251,7 @@ class RandomSetter
HashMapType* m_hashmaps;
SparseMatrixType* mp_target;
int m_outerPackets;
+ unsigned char m_keyBitsOffset;
};
#endif // EIGEN_RANDOMSETTER_H
diff --git a/Eigen/src/Sparse/SparseMatrixBase.h b/Eigen/src/Sparse/SparseMatrixBase.h
index b03fb5cee..2695f23b1 100644
--- a/Eigen/src/Sparse/SparseMatrixBase.h
+++ b/Eigen/src/Sparse/SparseMatrixBase.h
@@ -64,10 +64,11 @@ class SparseMatrixBase : public MatrixBase<Derived>
{
// std::cout << "Derived& operator=(const MatrixBase<OtherDerived>& other)\n";
const bool transpose = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit);
-// std::cout << "eval transpose = " << transpose << "\n";
+ ei_assert((!transpose) && "the transpose operation is supposed to be handled in SparseMatrix::operator=");
const int outerSize = other.outerSize();
- typedef typename ei_meta_if<transpose, LinkedVectorMatrix<Scalar,Flags&RowMajorBit>, Derived>::ret TempType;
- TempType temp(other.rows(), other.cols());
+ //typedef typename ei_meta_if<transpose, LinkedVectorMatrix<Scalar,Flags&RowMajorBit>, Derived>::ret TempType;
+ // thanks to shallow copies, we always eval to a tempary
+ Derived temp(other.rows(), other.cols());
temp.startFill(std::max(this->rows(),this->cols())*2);
for (int j=0; j<outerSize; ++j)
diff --git a/Eigen/src/Sparse/SparseSetter.h b/Eigen/src/Sparse/SparseSetter.h
deleted file mode 100644
index 20569d3ba..000000000
--- a/Eigen/src/Sparse/SparseSetter.h
+++ /dev/null
@@ -1,138 +0,0 @@
-// This file is part of Eigen, a lightweight C++ template library
-// for linear algebra. Eigen itself is part of the KDE project.
-//
-// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
-//
-// Eigen is free software; you can redistribute it and/or
-// modify it under the terms of the GNU Lesser General Public
-// License as published by the Free Software Foundation; either
-// version 3 of the License, or (at your option) any later version.
-//
-// Alternatively, you can redistribute it and/or
-// modify it under the terms of the GNU General Public License as
-// published by the Free Software Foundation; either version 2 of
-// the License, or (at your option) any later version.
-//
-// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
-// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
-// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
-// GNU General Public License for more details.
-//
-// You should have received a copy of the GNU Lesser General Public
-// License and a copy of the GNU General Public License along with
-// Eigen. If not, see <http://www.gnu.org/licenses/>.
-
-#ifndef EIGEN_SPARSESETTER_H
-#define EIGEN_SPARSESETTER_H
-
-template<typename MatrixType, int AccessPattern,
- int IsSupported = ei_support_access_pattern<MatrixType,AccessPattern>::ret>
-struct ei_sparse_setter_selector;
-
-/** \class SparseSetter
- *
- * Goal: provides a unified API to fill/update a dense or sparse matrix.
- *
- * Usage:
- * \code
- * {
- * SparseSetter<MatrixType, RandomAccessPattern> w(m);
- * for (...) w->coeffRef(rand(),rand()) = rand();
- * }
- * \endcode
- *
- * In the above example we want to fill a matrix m (could be a SparseMatrix or whatever other matrix type)
- * in a random fashion (whence the RandomAccessPattern). Internally, if \a MatrixType supports random writes
- * then \c w behaves as a pointer to m, and m is filled directly. Otherwise, a temporary matrix supporting
- * random writes is created and \c w behaves as a pointer to this temporary object. When the object \c w
- * is deleted (at the end of the block), then the temporary object is assigned to the matrix m.
- *
- * So far we can distinghished 4 types of access pattern:
- * - FullyCoherentAccessPattern (if col major, i+j*rows must increase)
- * - InnerCoherentAccessPattern (if col major, i must increase for each column j)
- * - OuterCoherentAccessPattern (if col major, the column j is set in a random order, but j must increase)
- * - RandomAccessPattern
- *
- * See the wiki for more details.
- *
- * The template class ei_support_access_pattern is used to determine the type of the temporary object (which
- * can be a reference to \a MatrixType if \a MatrixType support \a AccessPattern)
- *
- * Currently only the RandomAccessPattern seems to work as expected.
- *
- * \todo define the API for each kind of access pattern
- * \todo allows both update and set modes (set start a new matrix)
- * \todo implement the OuterCoherentAccessPattern
- *
- */
-template<typename MatrixType,
- int AccessPattern,
- typename WrapperType = typename ei_sparse_setter_selector<MatrixType,AccessPattern>::type>
-class SparseSetter
-{
- typedef typename ei_unref<WrapperType>::type _WrapperType;
- public:
-
- inline SparseSetter(MatrixType& matrix) : m_wrapper(matrix), mp_matrix(&matrix) {}
-
- ~SparseSetter()
- { *mp_matrix = m_wrapper; }
-
- inline _WrapperType* operator->() { return &m_wrapper; }
-
- inline _WrapperType& operator*() { return m_wrapper; }
-
- protected:
-
- WrapperType m_wrapper;
- MatrixType* mp_matrix;
-};
-
-template<typename MatrixType, int AccessPattern>
-struct ei_sparse_setter_selector<MatrixType, AccessPattern, AccessPatternSupported>
-{
- typedef MatrixType& type;
-};
-
-// forward each derived of SparseMatrixBase to the generic SparseMatrixBase specializations
-template<typename Scalar, int Flags, int AccessPattern>
-struct ei_sparse_setter_selector<SparseMatrix<Scalar,Flags>, AccessPattern, AccessPatternNotSupported>
-: public ei_sparse_setter_selector<SparseMatrixBase<SparseMatrix<Scalar,Flags> >,AccessPattern, AccessPatternNotSupported>
-{};
-
-template<typename Scalar, int Flags, int AccessPattern>
-struct ei_sparse_setter_selector<LinkedVectorMatrix<Scalar,Flags>, AccessPattern, AccessPatternNotSupported>
-: public ei_sparse_setter_selector<LinkedVectorMatrix<SparseMatrix<Scalar,Flags> >,AccessPattern, AccessPatternNotSupported>
-{};
-
-template<typename Scalar, int Flags, int AccessPattern>
-struct ei_sparse_setter_selector<HashMatrix<Scalar,Flags>, AccessPattern, AccessPatternNotSupported>
-: public ei_sparse_setter_selector<HashMatrix<SparseMatrix<Scalar,Flags> >,AccessPattern, AccessPatternNotSupported>
-{};
-
-// generic SparseMatrixBase specializations
-template<typename Derived>
-struct ei_sparse_setter_selector<SparseMatrixBase<Derived>, RandomAccessPattern, AccessPatternNotSupported>
-{
- typedef HashMatrix<typename Derived::Scalar, Derived::Flags> type;
-};
-
-template<typename Derived>
-struct ei_sparse_setter_selector<SparseMatrixBase<Derived>, OuterCoherentAccessPattern, AccessPatternNotSupported>
-{
- typedef HashMatrix<typename Derived::Scalar, Derived::Flags> type;
-};
-
-template<typename Derived>
-struct ei_sparse_setter_selector<SparseMatrixBase<Derived>, InnerCoherentAccessPattern, AccessPatternNotSupported>
-{
- typedef LinkedVectorMatrix<typename Derived::Scalar, Derived::Flags> type;
-};
-
-template<typename Derived>
-struct ei_sparse_setter_selector<SparseMatrixBase<Derived>, FullyCoherentAccessPattern, AccessPatternNotSupported>
-{
- typedef SparseMatrix<typename Derived::Scalar, Derived::Flags> type;
-};
-
-#endif // EIGEN_SPARSESETTER_H