const unsigned int dst_r=0,
const unsigned int dst_c=0) const;
- /**
+ /**
+ * Copy a subset of the rows and columns of another matrix into the
+ * current object.
+ *
+ * @param matrix The matrix from which a subset is to be taken from.
+ * @param row_index_set The set of rows of @p matrix from which to extract.
+ * @param column_index_set The set of columns of @p matrix from which to extract.
+ * @precondition The number of elements in @p row_index_set and
+ * @p column_index_set shall be equal to the number of
+ * rows and columns in the current object. In other words,
+ * the current object is not resized for this operation.
+ */
+ template <typename MatrixType>
+ void extract_submatrix_from (const MatrixType &matrix,
+ const std::vector<unsigned int> &row_index_set,
+ const std::vector<unsigned int> &column_index_set);
+
+ /**
+ * Copy the elements of the current matrix object into a specified
+ * set of rows and columns of another matrix. This thus a scatter operation.
+ *
+ * @param row_index_set The rows of @matrix into which to write.
+ * @param column_index_set The columns of @matrix into which to write.
+ * @param matrix The matrix within which certain elements are to be replaced.
+ * @precondition The number of elements in @p row_index_set and
+ * @p column_index_set shall be equal to the number of
+ * rows and columns in the current object. In other words,
+ * the current object is not resized for this operation.
+ */
+ template <typename MatrixType>
+ void
+ scatter_matrix_to (const std::vector<unsigned int> &row_index_set,
+ const std::vector<unsigned int> &column_index_set,
+ MatrixType &matrix) const;
+
+ /**
* Fill rectangular block.
*
* A rectangular block of the
const std::vector<unsigned int> &p_rows,
const std::vector<unsigned int> &p_cols);
+ /**
+ * Set a particular entry of the matrix to a value. Thus, calling
+ * <code>A.set(1,2,3.141);</code> is entirely equivalent to the operation
+ * <code>A(1,2) = 3.141;</code>. This function exists for compatibility
+ * with the various sparse matrix objects.
+ *
+ * @param i The row index of the element to be set.
+ * @param j The columns index of the element to be set.
+ * @param value The value to be written into the element.
+ */
+ void set (const unsigned int i,
+ const unsigned int j,
+ const number value);
/**
* @}
*/
}
+
template <typename number>
template <class MATRIX>
void
}
+
+template <typename number>
+template <typename MatrixType>
+inline
+void
+FullMatrix<number>::extract_submatrix_from (const MatrixType &matrix,
+ const std::vector<unsigned int> &row_index_set,
+ const std::vector<unsigned int> &column_index_set)
+{
+ AssertDimension(row_index_set.size(), this->n_rows());
+ AssertDimension(column_index_set.size(), this->n_cols());
+
+ const unsigned int n_rows_submatrix = row_index_set.size();
+ const unsigned int n_cols_submatrix = column_index_set.size();
+
+ for (unsigned int sub_row = 0; sub_row < n_rows_submatrix; ++sub_row)
+ for (unsigned int sub_col = 0; sub_col < n_cols_submatrix; ++sub_col)
+ (*this)(sub_row, sub_col) = matrix.el(row_index_set[sub_row], column_index_set[sub_col]);
+}
+
+
+
+template <typename number>
+template <typename MatrixType>
+inline
+void
+FullMatrix<number>::scatter_matrix_to (const std::vector<unsigned int> &row_index_set,
+ const std::vector<unsigned int> &column_index_set,
+ MatrixType &matrix) const
+{
+ AssertDimension(row_index_set.size(), this->n_rows());
+ AssertDimension(column_index_set.size(), this->n_cols());
+
+ const unsigned int n_rows_submatrix = row_index_set.size();
+ const unsigned int n_cols_submatrix = column_index_set.size();
+
+ for (unsigned int sub_row = 0; sub_row < n_rows_submatrix; ++sub_row)
+ for (unsigned int sub_col = 0; sub_col < n_cols_submatrix; ++sub_col)
+ matrix.set(row_index_set[sub_row],
+ column_index_set[sub_col],
+ (*this)(sub_row, sub_col));
+}
+
+
+template <typename number>
+inline
+void
+FullMatrix<number>::set (const unsigned int i,
+ const unsigned int j,
+ const number value)
+{
+ (*this)(i,j) = value;
+}
+
+
+
template <typename number>
template<typename number2>
void
--- /dev/null
+//---------------------------- full_matrix_05.cc,v ---------------------------
+// $Id$
+//
+// Copyright (C) 2011, 2012 by the deal.II authors
+//
+// This file is subject to QPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//---------------------------- full_matrix_05.cc,v ---------------------------
+
+// check method FullMatrix::extract_submatrix_from
+
+#include "../tests.h"
+#include <cmath>
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <cstdlib>
+
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/full_matrix.h>
+
+void test ()
+{
+ // create a matrix with known
+ // elements
+ FullMatrix<double> A(10,12);
+ for (unsigned int i=0; i<A.m(); ++i)
+ for (unsigned int j=0; j<A.n(); ++j)
+ A(i,j) = i+j;
+
+ // pick every other row and column
+ std::vector<unsigned int> rows (A.m()/2);
+ for (unsigned int i=0; i<rows.size(); ++i)
+ rows[i] = 2*i;
+
+ std::vector<unsigned int> cols (A.n()/2);
+ for (unsigned int i=0; i<cols.size(); ++i)
+ cols[i] = 2*i;
+
+ // do the extraction
+ FullMatrix<double> X(rows.size(), cols.size());
+ X.extract_submatrix_from (A, rows, cols);
+
+ // verify that the elements are
+ // correct
+ for (unsigned int i=0; i<X.m(); ++i)
+ for (unsigned int j=0; j<X.n(); ++j)
+ Assert (X(i,j) == 2*i + 2*j,
+ ExcInternalError());
+
+ deallog << "OK" << std::endl;
+}
+
+
+int
+main ()
+{
+ const std::string logname = JobIdentifier::base_name(__FILE__) + std::string("/output");
+ std::ofstream logfile(logname.c_str());
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ test();
+}
--- /dev/null
+//---------------------------- full_matrix_06.cc,v ---------------------------
+// $Id$
+//
+// Copyright (C) 2011, 2012 by the deal.II authors
+//
+// This file is subject to QPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//---------------------------- full_matrix_06.cc,v ---------------------------
+
+// check method FullMatrix::scatter_matrix_to
+
+#include "../tests.h"
+#include <cmath>
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <cstdlib>
+
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/full_matrix.h>
+
+void test ()
+{
+ // create a matrix with known
+ // elements
+ FullMatrix<double> A(5,6);
+ for (unsigned int i=0; i<A.m(); ++i)
+ for (unsigned int j=0; j<A.n(); ++j)
+ A(i,j) = i+j;
+
+ // pick every other row and column
+ std::vector<unsigned int> rows (A.m());
+ for (unsigned int i=0; i<rows.size(); ++i)
+ rows[i] = 2*i;
+
+ std::vector<unsigned int> cols (A.n());
+ for (unsigned int i=0; i<cols.size(); ++i)
+ cols[i] = 2*i;
+
+ // do the scatter
+ FullMatrix<double> X(rows.size()*2, cols.size()*2);
+ A.scatter_matrix_to (rows, cols, X);
+
+ // verify that the elements are
+ // correct
+ for (unsigned int i=0; i<X.m(); ++i)
+ for (unsigned int j=0; j<X.n(); ++j)
+ if ((i % 2 == 0) && (j % 2 == 0))
+ Assert (X(i,j) == i/2 + j/2,
+ ExcInternalError())
+ else
+ Assert (X(i,j) == 0,
+ ExcInternalError());
+
+ deallog << "OK" << std::endl;
+}
+
+
+int
+main ()
+{
+ const std::string logname = JobIdentifier::base_name(__FILE__) + std::string("/output");
+ std::ofstream logfile(logname.c_str());
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ test();
+}