<h3>lac</h3>
<ol>
+ <li> <p> New: Functions <code
+ class="member">SparsityPattern::copy_from</code> and <code
+ class="member">SparseMatrix::copy_from</code> allow to copy a
+ full matrix into a sparse matrix.
+ <br>
+ (WB 2002/02/06)
+ </p>
</ol>
template SparseMatrix<TYPEMAT> &
SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+template
+void SparseMatrix<TYPEMAT>::copy_from (const FullMatrix<TYPE2> &);
+
template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
const SparseMatrix<TYPE2> &);
#include <lac/sparsity_pattern.h>
template<typename number> class Vector;
+template<typename number> class FullMatrix;
/**
* Sparse matrix.
template <typename ForwardIterator>
void copy_from (const ForwardIterator begin,
const ForwardIterator end);
+
+ /**
+ * Copy the nonzero entries of a
+ * full matrix into this
+ * object. Previous content is
+ * deleted. Note that the
+ * underlying sparsity pattern
+ * must be appropriate to hold
+ * the nonzero entries of the
+ * full matrix.
+ */
+ template <typename somenumber>
+ void copy_from (const FullMatrix<somenumber> &matrix);
/**
* Add @p{matrix} scaled by
#include <lac/sparse_matrix.h>
#include <lac/vector.h>
+#include <lac/full_matrix.h>
// we only need output streams, but older compilers did not provide
};
+
+template <typename number>
+template <typename somenumber>
+void
+SparseMatrix<number>::copy_from (const FullMatrix<somenumber> &matrix)
+{
+ // first delete previous content
+ clear ();
+
+ // then copy old matrix
+ for (unsigned int row=0; row<matrix.m(); ++row)
+ for (unsigned int col=0; col<matrix.n(); ++col)
+ if (matrix(row,col) != 0)
+ set (row, col, matrix(row,col));
+};
+
+
+
template <typename number>
template <typename somenumber>
void
#include <vector>
+template <typename number> class FullMatrix;
template <typename number> class SparseMatrix;
class CompressedSparsityPattern;
/**
* Copy data from an object of
* type
- * @ref{CompressedSparsityPattern}. Previous
- * content of this object is
- * lost.
+ * @ref{CompressedSparsityPattern}.
+ * Previous content of this
+ * object is lost, and the
+ * sparsity pattern is in
+ * compressed mode afterwards.
*/
void copy_from (const CompressedSparsityPattern &csp);
+
+ /**
+ * Take a full matrix and use its
+ * nonzero entries to generate a
+ * sparse matrix entry pattern
+ * for this object.
+ *
+ * Previous content of this
+ * object is lost, and the
+ * sparsity pattern is in
+ * compressed mode afterwards.
+ */
+ template <typename number>
+ void copy_from (const FullMatrix<number> &matrix);
/**
* Return whether the object is empty. It
#include <lac/sparsity_pattern.h>
+#include <lac/full_matrix.h>
#include <lac/compressed_sparsity_pattern.h>
#include <iostream>
};
+template <typename number>
+void SparsityPattern::copy_from (const FullMatrix<number> &matrix)
+{
+ // first init with the number of
+ // entries per row
+ std::vector<unsigned int> entries_per_row (matrix.m(), 0);
+ for (unsigned int row=0; row<matrix.m(); ++row)
+ for (unsigned int col=0; col<matrix.n(); ++col)
+ if (matrix(row,col) != 0)
+ ++entries_per_row[row];
+ reinit (matrix.m(), matrix.n(), entries_per_row);
+
+ // now set entries
+ for (unsigned int row=0; row<matrix.m(); ++row)
+ for (unsigned int col=0; col<matrix.n(); ++col)
+ if (matrix(row,col) != 0)
+ add (row,col);
+
+ // finally compress
+ compress ();
+};
+
+
bool
SparsityPattern::empty () const
max_dim * sizeof(unsigned int) +
max_vec_len * sizeof(unsigned int));
};
+
+
+
+// explicit instantiations
+template void SparsityPattern::copy_from (const FullMatrix<float> &);
+template void SparsityPattern::copy_from (const FullMatrix<double> &);