]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add identity matrix class.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 24 Apr 2006 20:24:27 +0000 (20:24 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 24 Apr 2006 20:24:27 +0000 (20:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@12867 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/identity_matrix.h [new file with mode: 0644]

diff --git a/deal.II/lac/include/lac/identity_matrix.h b/deal.II/lac/include/lac/identity_matrix.h
new file mode 100644 (file)
index 0000000..a244633
--- /dev/null
@@ -0,0 +1,275 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2006 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.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__identity_matrix_h
+#define __deal2__identity_matrix_h
+
+
+#include <base/config.h>
+#include <lac/exceptions.h>
+
+
+/*! @addtogroup Matrix1
+ *@{
+ */
+
+
+/**
+ * Implementation of a simple class representing the identity matrix
+ * of a given size, i.e. a matrix with entries
+ * $A_{ij}=\delta_{ij}$. While it has the most important ingredients
+ * of a matrix, in particular that one can ask for its size and
+ * perform matrix-vector products with it, a matrix of this type is
+ * really only useful in two contexts: preconditioning and
+ * initializing other matrices.
+ *
+ * <h4>Initialization</h4>
+ *
+ * The main usefulness of this class lies in its ability to initialize
+ * other matrix, like this:
+ @verbatim
+   FullMatrix<double> identity (IdentityMatrix(10));
+ @endverbatim
+ *
+ * This creates a $10\times 10$ matrix with ones on the diagonal and
+ * zeros everywhere else. All matrix types have conversion
+ * constructors and can therefore be filled rather easily with
+ * identity matrices.
+ *
+ *
+ * <h4>Preconditioning</h4>
+ *
+ * No preconditioning at all is equivalent to preconditioning with
+ * preconditioning with the identity matrix. deal.II has a specialized
+ * class for this purpose, PreconditionIdentity, than can be used in a
+ * context as shown in the documentation of that class. The present
+ * class can be used in much the same way, although without any
+ * additional benefit:
+ @verbatim
+  SolverControl           solver_control (1000, 1e-12);
+  SolverCG<>              cg (solver_control);
+  cg.solve (system_matrix, solution, system_rhs,
+           IdentityMatrix(solution.size()));
+ @endverbatim 
+ *
+ *
+ * @author Wolfgang Bangerth, 2006
+ */
+class IdentityMatrix
+{
+  public:
+                                    /**
+                                     * Default constructor. Creates a
+                                     * zero-sized matrix that should
+                                     * be resized later on using the
+                                     * reinit() function.
+                                     */
+    IdentityMatrix ();
+
+                                    /**
+                                     * Constructor. Creates a
+                                     * identity matrix of size #n.
+                                     */
+    IdentityMatrix (const unsigned int n);
+
+                                    /**
+                                     * Resize the matrix to be of
+                                     * size #n by #n.
+                                     */
+    void reinit (const unsigned int n);
+
+                                    /**
+                                     * Number of rows of this
+                                     * matrix. For the present
+                                     * matrix, the number of rows and
+                                     * columns are equal, of course.
+                                     */
+    unsigned int m () const;
+    
+                                    /**
+                                     * Number of columns of this
+                                     * matrix. For the present
+                                     * matrix, the number of rows and
+                                     * columns are equal, of course.
+                                     */
+    unsigned int n () const;
+
+                                    /**
+                                     * Matrix-vector
+                                     * multiplication. For the
+                                     * present case, this of course
+                                     * amounts to simply copying the
+                                     * input vector to the output
+                                     * vector.
+                                     */
+    template <class VECTOR1, class VECTOR2>
+    void vmult (VECTOR1       &out,
+               const VECTOR2 &in) const;
+
+                                    /**
+                                     * Matrix-vector multiplication
+                                     * with addition to the output
+                                     * vector. For the present case,
+                                     * this of course amounts to
+                                     * simply adding the input
+                                     * vector to the output vector.
+                                     */
+    template <class VECTOR1, class VECTOR2>
+    void vmult_add (VECTOR1       &out,
+                   const VECTOR2 &in) const;
+    
+                                    /**
+                                     * Matrix-vector multiplication
+                                     * with the transpose matrix. For
+                                     * the present case, this of
+                                     * course amounts to simply
+                                     * copying the input vector to
+                                     * the output vector.
+                                     */
+    template <class VECTOR1, class VECTOR2>
+    void Tvmult (VECTOR1       &out,
+                const VECTOR2 &in) const;
+    
+
+                                    /**
+                                     * Matrix-vector multiplication
+                                     * with the transpose matrix,
+                                     * with addition to the output
+                                     * vector. For the present case,
+                                     * this of course amounts to
+                                     * simply adding the input vector
+                                     * to the output vector.
+                                     */
+    template <class VECTOR1, class VECTOR2>
+    void Tvmult_add (VECTOR1       &out,
+                    const VECTOR2 &in) const;
+  private:
+    
+                                    /**
+                                     * Number of rows and columns of
+                                     * this matrix.
+                                     */
+    unsigned int size;
+};
+
+
+
+
+// ------------------------- inline and template functions -------------
+#ifndef DOXYGEN
+
+
+inline
+IdentityMatrix::IdentityMatrix ()
+               :
+               size (0)
+{}
+
+
+
+inline
+IdentityMatrix::IdentityMatrix (const unsigned int n)
+               :
+               size (n)
+{}
+
+
+
+inline
+void
+IdentityMatrix::reinit (const unsigned int n)
+{
+  size = n;
+}
+
+
+
+inline
+unsigned int
+IdentityMatrix::m () const
+{
+  return size;
+}
+
+
+
+inline
+unsigned int
+IdentityMatrix::n () const
+{
+  return size;
+}
+
+
+
+template <class VECTOR1, class VECTOR2>
+inline
+void
+IdentityMatrix::vmult (VECTOR1       &out,
+                      const VECTOR2 &in) const
+{
+  Assert (out.size() == size, ExcDimensionMismatch (out.size(), size));
+  Assert (in.size() == size, ExcDimensionMismatch (in.size(), size));
+
+  out = in;
+}
+
+
+
+template <class VECTOR1, class VECTOR2>
+inline
+void
+IdentityMatrix::vmult_add (VECTOR1       &out,
+                          const VECTOR2 &in) const
+{
+  Assert (out.size() == size, ExcDimensionMismatch (out.size(), size));
+  Assert (in.size() == size, ExcDimensionMismatch (in.size(), size));
+
+  out += in;
+}
+
+
+
+template <class VECTOR1, class VECTOR2>
+inline
+void
+IdentityMatrix::Tvmult (VECTOR1       &out,
+                       const VECTOR2 &in) const
+{
+  Assert (out.size() == size, ExcDimensionMismatch (out.size(), size));
+  Assert (in.size() == size, ExcDimensionMismatch (in.size(), size));
+
+  out = in;
+}
+
+
+
+template <class VECTOR1, class VECTOR2>
+inline
+void
+IdentityMatrix::Tvmult_add (VECTOR1       &out,
+                           const VECTOR2 &in) const
+{
+  Assert (out.size() == size, ExcDimensionMismatch (out.size(), size));
+  Assert (in.size() == size, ExcDimensionMismatch (in.size(), size));
+
+  out += in;
+}
+
+
+#endif
+
+/**@}*/
+
+#endif
+

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.