From: bangerth Date: Tue, 24 Apr 2012 09:24:25 +0000 (+0000) Subject: New functions in FullMatrix: extract_submatrix_from, scatter_matrix_to, set. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=275f72f4369ed55467974e0448cf340d1dde1a21;p=dealii-svn.git New functions in FullMatrix: extract_submatrix_from, scatter_matrix_to, set. git-svn-id: https://svn.dealii.org/trunk@25434 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 7aa1778a3b..fbbea4a191 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -165,6 +165,11 @@ enabled due to a missing include file in file

Specific improvements

    +
  1. New: The FullMatrix::extract_submatrix_from, FullMatrix::scatter_matrix_to, +FullMatrix::set functions are new. +
    +(Andrew McBride, Jean Paul Pelteret, Wolfgang Bangerth, 2012/04/24) +
  2. Fixed: The method FEValues::inverse_jacobian() previously returned the transpose of the inverse Jacobians instead of just the inverse Jacobian as documented. This is now fixed. diff --git a/deal.II/include/deal.II/lac/full_matrix.h b/deal.II/include/deal.II/lac/full_matrix.h index 7c049d7d87..cde6d5b188 100644 --- a/deal.II/include/deal.II/lac/full_matrix.h +++ b/deal.II/include/deal.II/lac/full_matrix.h @@ -421,7 +421,42 @@ class FullMatrix : public Table<2,number> 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 + void extract_submatrix_from (const MatrixType &matrix, + const std::vector &row_index_set, + const std::vector &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 + void + scatter_matrix_to (const std::vector &row_index_set, + const std::vector &column_index_set, + MatrixType &matrix) const; + + /** * Fill rectangular block. * * A rectangular block of the @@ -477,6 +512,19 @@ class FullMatrix : public Table<2,number> const std::vector &p_rows, const std::vector &p_cols); + /** + * Set a particular entry of the matrix to a value. Thus, calling + * A.set(1,2,3.141); is entirely equivalent to the operation + * A(1,2) = 3.141;. 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); /** * @} */ @@ -1453,6 +1501,7 @@ FullMatrix::copy_from (const MATRIX& M) } + template template void @@ -1466,6 +1515,62 @@ FullMatrix::copy_transposed (const MATRIX& M) } + +template +template +inline +void +FullMatrix::extract_submatrix_from (const MatrixType &matrix, + const std::vector &row_index_set, + const std::vector &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 +template +inline +void +FullMatrix::scatter_matrix_to (const std::vector &row_index_set, + const std::vector &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 +inline +void +FullMatrix::set (const unsigned int i, + const unsigned int j, + const number value) +{ + (*this)(i,j) = value; +} + + + template template void diff --git a/tests/lac/full_matrix_05.cc b/tests/lac/full_matrix_05.cc new file mode 100644 index 0000000000..5685759a48 --- /dev/null +++ b/tests/lac/full_matrix_05.cc @@ -0,0 +1,67 @@ +//---------------------------- 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 +#include +#include +#include +#include + +#include +#include + +void test () +{ + // create a matrix with known + // elements + FullMatrix A(10,12); + for (unsigned int i=0; i rows (A.m()/2); + for (unsigned int i=0; i cols (A.n()/2); + for (unsigned int i=0; i X(rows.size(), cols.size()); + X.extract_submatrix_from (A, rows, cols); + + // verify that the elements are + // correct + for (unsigned int i=0; i +#include +#include +#include +#include + +#include +#include + +void test () +{ + // create a matrix with known + // elements + FullMatrix A(5,6); + for (unsigned int i=0; i rows (A.m()); + for (unsigned int i=0; i cols (A.n()); + for (unsigned int i=0; i 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