]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New functions in FullMatrix: extract_submatrix_from, scatter_matrix_to, set.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 24 Apr 2012 09:24:25 +0000 (09:24 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 24 Apr 2012 09:24:25 +0000 (09:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@25434 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/lac/full_matrix.h
tests/lac/full_matrix_05.cc [new file with mode: 0644]
tests/lac/full_matrix_05/cmp/generic [new file with mode: 0644]
tests/lac/full_matrix_06.cc [new file with mode: 0644]
tests/lac/full_matrix_06/cmp/generic [new file with mode: 0644]

index 7aa1778a3b02bb652253639305ef336028894620..fbbea4a1919888567f83d8da4af17d4a375fb24d 100644 (file)
@@ -165,6 +165,11 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: The FullMatrix::extract_submatrix_from, FullMatrix::scatter_matrix_to,
+FullMatrix::set functions are new.
+<br>
+(Andrew McBride, Jean Paul Pelteret, Wolfgang Bangerth, 2012/04/24)
+
 <li> Fixed:
 The method FEValues<dim>::inverse_jacobian() previously returned the transpose of the
 inverse Jacobians instead of just the inverse Jacobian as documented. This is now fixed.
index 7c049d7d879b024a57aa9fbafa675700a71fd94c..cde6d5b188709cd7ffcd815a4c5c8837bb36545a 100644 (file)
@@ -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 <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
@@ -477,6 +512,19 @@ class FullMatrix : public Table<2,number>
                            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);
     /**
      * @}
      */
@@ -1453,6 +1501,7 @@ FullMatrix<number>::copy_from (const MATRIX& M)
 }
 
 
+
 template <typename number>
 template <class MATRIX>
 void
@@ -1466,6 +1515,62 @@ FullMatrix<number>::copy_transposed (const MATRIX& M)
 }
 
 
+
+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
diff --git a/tests/lac/full_matrix_05.cc b/tests/lac/full_matrix_05.cc
new file mode 100644 (file)
index 0000000..5685759
--- /dev/null
@@ -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 <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();
+}
diff --git a/tests/lac/full_matrix_05/cmp/generic b/tests/lac/full_matrix_05/cmp/generic
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK
diff --git a/tests/lac/full_matrix_06.cc b/tests/lac/full_matrix_06.cc
new file mode 100644 (file)
index 0000000..e5e1978
--- /dev/null
@@ -0,0 +1,71 @@
+//----------------------------  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();
+}
diff --git a/tests/lac/full_matrix_06/cmp/generic b/tests/lac/full_matrix_06/cmp/generic
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.