]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new class PoiterMatrixAux
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 2 Feb 2006 09:29:11 +0000 (09:29 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 2 Feb 2006 09:29:11 +0000 (09:29 +0000)
git-svn-id: https://svn.dealii.org/trunk@12223 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.html
deal.II/lac/include/lac/block_matrix_array.h
deal.II/lac/include/lac/pointer_matrix.h

index c14604d344b8fb0458073e0894a30d2d2aaf3da8..0576cb49bfa21041e42f64d10265715a792f9ab5 100644 (file)
@@ -298,6 +298,23 @@ inconvenience this causes.
 <h3>lac</h3>
 
 <ol>
+
+
+  >li> <p> Improved: <code class="class">BlockMatrixArray</class>::<code
+  class="member">enter_aux</class> allows using matrices without adding vector
+  muiltplication by using <code class="class">PointerMatrixAux</class>.
+  <br>
+  (GK 2006/02/02)
+  </p>  
+
+  <li> <p> New: The class <code class="class">PointerMatrixAux</class> was
+  introduced for use with matrix classes lacking the adding vector
+  multiplication functions. It implements these by using its own auxiliary
+  vector.
+  <br>
+  (GK 2006/02/02)
+  </p>
+
   <li> <p>
        Improved: The <code class="class">FilteredMatrix</code> class
        was able to filter only the <code class="class">SparseMatrix</code>
index 4a867ecfd533d63d215b7e2238a318a965ce66cc..87320b854b9342731494a6d6dbce6330c5510b3b 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -185,6 +185,21 @@ class BlockMatrixArray : public Subscriptor
                const unsigned int col,
                const double       prefix = 1.,
                const bool         transpose = false);
+
+                                    /**
+                                     * Add an entry like with enter,
+                                     * but use PointerMatrixAux for
+                                     * matrices not having functions
+                                     * vmult_add() and TVmult_add().
+                                     */
+    template <class MATRIX>
+    void enter_aux (VectorMemory<Vector<number> >& mem,
+                   const MATRIX      &matrix,
+                   const unsigned int row,
+                   const unsigned int col,
+                   const double       prefix = 1.,
+                   const bool         transpose = false);
+
     
                                     /**
                                      * Delete all entries, i.e. reset
@@ -284,13 +299,30 @@ class BlockMatrixArray : public Subscriptor
     {
       public:
                                         /**
-                                         * Constructor initializing all
-                                         * data fields.
+                                         * Constructor initializing
+                                         * all data fields. A
+                                         * PointerMatrix object is
+                                         * generated for
+                                         * <tt>matrix</tt>.
                                          */
        template<class MATRIX> 
        Entry (const MATRIX& matrix,
               unsigned row, unsigned int col,
               double prefix, bool transpose);
+       
+                                        /**
+                                         * Constructor initializing
+                                         * all data fields. A
+                                         * PointerMatrixAux object is
+                                         * generated for
+                                         * <tt>matrix</tt>.
+                                         */
+       template<class MATRIX> 
+       Entry (const MATRIX& matrix,
+              unsigned row, unsigned int col,
+              double prefix, bool transpose,
+              VectorMemory<Vector<number> >& mem);
+       
                                         /**
                                          * Copy constructor
                                          * invalidating the old
@@ -484,6 +516,21 @@ class BlockTrianglePrecondition
                const double       prefix = 1.,
                const bool         transpose = false);
 
+                                    /**
+                                     * Enter a block. This calls
+                                     * BlockMatrixArray::enter_aux(). Remember
+                                     * that the diagonal blocks
+                                     * should actually be inverse
+                                     * matrices or preconditioners.
+                                     */
+    template <class MATRIX>
+    void enter_aux (VectorMemory<Vector<double> >& mem,
+                   const MATRIX      &matrix,
+                   const unsigned int row,
+                   const unsigned int col,
+                   const double       prefix = 1.,
+                   const bool         transpose = false);
+
                                      /**
                                      * Preconditioning.
                                      */
@@ -578,15 +625,38 @@ class BlockTrianglePrecondition
 template <typename number>
 template <class MATRIX>
 inline
-BlockMatrixArray<number>::Entry::Entry (const MATRIX& matrix,
-                                       unsigned row, unsigned int col,
-                                       double prefix, bool transpose)
+BlockMatrixArray<number>::Entry::Entry (
+  const MATRIX& m,
+  unsigned int row,
+  unsigned int col,
+  double prefix,
+  bool transpose)
                :
                row (row),
                col (col),
                prefix (prefix),
                transpose (transpose),
-               matrix (new PointerMatrix<MATRIX, Vector<number> >(&matrix, typeid(*this).name()))
+               matrix (new PointerMatrix<MATRIX, Vector<number> >(&m, typeid(*this).name()))
+{}
+
+
+
+template <typename number>
+template <class MATRIX>
+inline
+BlockMatrixArray<number>::Entry::Entry (
+  const MATRIX& m,
+  unsigned row,
+  unsigned int col,
+  double prefix,
+  bool transpose,
+  VectorMemory<Vector<number> >& mem)
+               :
+               row (row),
+               col (col),
+               prefix (prefix),
+               transpose (transpose),
+               matrix (new PointerMatrixAux<MATRIX, Vector<number> >(&mem, &m, typeid(*this).name()))
 {}
 
 
@@ -595,9 +665,12 @@ template <typename number>
 template <class MATRIX>
 inline
 void
-BlockMatrixArray<number>::enter (const MATRIX& matrix,
-                                unsigned row, unsigned int col,
-                                double prefix, bool transpose)
+BlockMatrixArray<number>::enter (
+  const MATRIX& matrix,
+  unsigned int row,
+  unsigned int col,
+  double prefix,
+  bool transpose)
 {
   Assert (mem != 0, ExcNotInitialized());
   Assert(row<n_block_rows(), ExcIndexRange(row, 0, n_block_rows()));
@@ -606,6 +679,24 @@ BlockMatrixArray<number>::enter (const MATRIX& matrix,
 }
 
 
+template <typename number>
+template <class MATRIX>
+inline
+void
+BlockMatrixArray<number>::enter_aux (
+  VectorMemory<Vector<number> >& mem,
+  const MATRIX& matrix,
+  unsigned int row,
+  unsigned int col,
+  double prefix,
+  bool transpose)
+{
+  Assert(row<n_block_rows(), ExcIndexRange(row, 0, n_block_rows()));
+  Assert(col<n_block_cols(), ExcIndexRange(col, 0, n_block_cols()));
+  entries.push_back(Entry(matrix, row, col, prefix, transpose, mem));
+}
+
+
 template<typename number>
 struct BlockMatrixArrayPointerMatrixLess
 {
@@ -704,6 +795,23 @@ BlockTrianglePrecondition<number>::enter (const MATRIX& matrix,
 
 
 
+template <typename number>
+template <class MATRIX>
+inline
+void
+BlockTrianglePrecondition<number>::enter_aux (
+  VectorMemory<Vector<double> >& mem,
+  const MATRIX& matrix,
+  unsigned int row,
+  unsigned int col,
+  double prefix,
+  bool transpose)
+{
+  BlockMatrixArray<number>::enter_aux(mem, matrix, row, col, prefix, transpose);
+}
+
+
+
 #endif // DOXYGEN
 
 #endif
index 690da5a9c56fcc892a2057461c518eabf3b07d59..c38ea22c5d4b722c9dc2dcf6e792803062add963 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2002, 2003, 2004, 2005 by the deal.II authors
+//    Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -16,6 +16,8 @@
 #include <base/subscriptor.h>
 #include <base/smartpointer.h>
 
+template<class VECTOR> class VectorMemory;
+
 /*! @addtogroup Matrix2
  *@{
  */
@@ -233,6 +235,122 @@ class PointerMatrix : public PointerMatrixBase<VECTOR>
 };
 
 
+/**
+ * A pointer to be used as a matrix.  This class stores a pointer to a
+ * matrix and can be used as a matrix itself in iterative methods.
+ *
+ * The main purpose for the existence of this class is its base class,
+ * which only has a vector as template argument. Therefore, this
+ * interface provides an abstract base class for matrices.
+ *
+ * This class differs form PointerMatrix by its additional
+ * VectorMemory object and by the fact that it implements the
+ * functions vmult_add() and Tvmult_add() only using vmult() and
+ * Tvmult() of the MATRIX.
+ *
+ * @author Guido Kanschat 2006
+ */
+template<class MATRIX, class VECTOR>
+class PointerMatrixAux : public PointerMatrixBase<VECTOR>
+{
+  public:
+                                    /**
+                                     * Constructor.  The pointer in the
+                                     * argument is stored in this
+                                     * class. As usual, the lifetime of
+                                     * <tt>*M</tt> must be longer than the
+                                     * one of the PointerMatrixAux.
+                                     *
+                                     * If <tt>M</tt> is zero, no
+                                     * matrix is stored.
+                                     */
+    PointerMatrixAux (VectorMemory<VECTOR>* mem,
+                     const MATRIX* M=0);
+
+                                    /**
+                                     * Constructor. The name argument
+                                     * is used to identify the
+                                     * SmartPointer for this object.
+                                     */
+    PointerMatrixAux(VectorMemory<VECTOR>* mem,
+                    const char* name);
+    
+                                    /**
+                                     * Constructor. <tt>M</tt> points
+                                     * to a matrix which must live
+                                     * longer than the
+                                     * PointerMatrixAux. The name
+                                     * argument is used to identify
+                                     * the SmartPointer for this
+                                     * object.
+                                     */
+    PointerMatrixAux(VectorMemory<VECTOR>* mem,
+                    const MATRIX* M,
+                    const char* name);
+
+                                    // Use doc from base class
+    virtual void clear();
+    
+                                    /**
+                                     * Return whether the object is
+                                     * empty. 
+                                     */
+    bool empty () const;
+    
+                                    /**
+                                     * Assign a new matrix
+                                     * pointer. Deletes the old pointer
+                                     * and releases its matrix.
+                                     * @see SmartPointer
+                                     */
+    const PointerMatrixAux& operator= (const MATRIX* M);
+    
+                                    /**
+                                     * Matrix-vector product.
+                                     */
+    virtual void vmult (VECTOR& dst,
+                       const VECTOR& src) const;
+    
+                                    /**
+                                     * Tranposed matrix-vector product.
+                                     */
+    virtual void Tvmult (VECTOR& dst,
+                        const VECTOR& src) const;
+    
+                                    /**
+                                     * Matrix-vector product, adding to
+                                     * <tt>dst</tt>.
+                                     */
+    virtual void vmult_add (VECTOR& dst,
+                           const VECTOR& src) const;
+    
+                                    /**
+                                     * Tranposed matrix-vector product,
+                                     * adding to <tt>dst</tt>.
+                                     */
+    virtual void Tvmult_add (VECTOR& dst,
+                            const VECTOR& src) const;
+    
+  private:
+                                    /**
+                                     * Return the address of the
+                                     * matrix for comparison.
+                                     */
+    virtual const void* get() const;
+
+                                    /**
+                                     * Object for getting the
+                                     * auxiliary vector.
+                                     */
+    SmartPointer<VectorMemory<VECTOR> > mem;
+    
+                                    /**
+                                     * The pointer to the actual matrix.
+                                     */
+    SmartPointer<const MATRIX> m;
+};
+
+
 /**
  * This function helps you creating a PointerMatrixBase object if you
  * do not want to provide the full template arguments of
@@ -418,4 +536,117 @@ PointerMatrix<MATRIX, VECTOR>::get () const
 }
 
 
+//----------------------------------------------------------------------//
+
+
+template<class MATRIX, class VECTOR>
+PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
+  VectorMemory<VECTOR>* mem,
+  const MATRIX* M)
+               : mem(mem, typeid(*this).name()),
+                 m(M, typeid(*this).name())
+{}
+
+
+template<class MATRIX, class VECTOR>
+PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
+  VectorMemory<VECTOR>* mem,
+  const char* name)
+               : mem(mem, name),
+                 m(0, name)
+{}
+
+
+template<class MATRIX, class VECTOR>
+PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
+  VectorMemory<VECTOR>* mem,
+  const MATRIX* M,
+  const char* name)
+               : mem(mem, name),
+                 m(M, name)
+{}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::clear ()
+{
+  m = 0;
+}
+
+
+template<class MATRIX, class VECTOR>
+inline const PointerMatrixAux<MATRIX, VECTOR>&
+PointerMatrixAux<MATRIX, VECTOR>::operator= (const MATRIX* M)
+{
+  m = M;
+  return *this;
+}
+
+
+template<class MATRIX, class VECTOR>
+inline bool
+PointerMatrixAux<MATRIX, VECTOR>::empty () const
+{
+  if (m == 0)
+    return true;
+  return m->empty();
+}
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::vmult (VECTOR& dst,
+                                     const VECTOR& src) const
+{
+  Assert (m != 0, ExcNotInitialized());
+  m->vmult (dst, src);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::Tvmult (VECTOR& dst,
+                                      const VECTOR& src) const
+{
+  Assert (m != 0, ExcNotInitialized());
+  m->Tvmult (dst, src);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::vmult_add (VECTOR& dst,
+                                         const VECTOR& src) const
+{
+  Assert (m != 0, ExcNotInitialized());
+  VECTOR* v = mem->alloc();
+  v->reinit(dst);
+  m->vmult (*v, src);
+  dst += *v;
+  mem->free(v);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::Tvmult_add (VECTOR& dst,
+                                             const VECTOR& src) const
+{
+  Assert (m != 0, ExcNotInitialized());
+  VECTOR* v = mem->alloc();
+  v->reinit(dst);
+  m->Tvmult (*v, src);
+  dst += *v;
+  mem->free(v);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline const void*
+PointerMatrixAux<MATRIX, VECTOR>::get () const
+{
+  return m;
+}
+
+
 #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.