]> https://gitweb.dealii.org/ - dealii.git/commitdiff
use new GrowingVectorMemory functionality
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 8 Jul 2010 23:41:10 +0000 (23:41 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 8 Jul 2010 23:41:10 +0000 (23:41 +0000)
git-svn-id: https://svn.dealii.org/trunk@21465 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/block_matrix_array.h
deal.II/lac/source/block_matrix_array.cc

index 93897ba88457d95eca2d81fe31eb9701cbf49712..5f0d55130aba93d01f73ed4121c43903cd4263e1 100644 (file)
@@ -1,6 +1,5 @@
 //---------------------------------------------------------------------------
 //    $Id$
-//    Version: $Name$
 //
 //    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
@@ -117,7 +116,8 @@ template <typename> class Vector;
  * BlockTrianglePrecondition.
  *
  * @see @ref GlossBlockLA "Block (linear algebra)"
- * @author Guido Kanschat, 2000-2005
+ * @author Guido Kanschat
+ * @date 2000-2005, 2010
  */
 template <typename number = double>
 class BlockMatrixArray : public Subscriptor
@@ -135,6 +135,29 @@ class BlockMatrixArray : public Subscriptor
                                      * Constructor fixing the
                                      * dimensions.
                                      */
+    BlockMatrixArray (const unsigned int n_block_rows,
+                     const unsigned int n_block_cols);
+
+                                    /**
+                                     * Initialize object
+                                     * completely. This is the
+                                     * function to call for an object
+                                     * created by the default
+                                     * constructor.
+                                     *
+                                     * @deprecated the last argument
+                                     * is ignored.
+                                     */
+    void initialize (const unsigned int n_block_rows,
+                    const unsigned int n_block_cols);
+    
+                                    /**
+                                     * Constructor fixing the
+                                     * dimensions.
+                                     *
+                                     * @deprecated the last argument
+                                     * is ignored.
+                                     */
     BlockMatrixArray (const unsigned int n_block_rows,
                      const unsigned int n_block_cols,
                      VectorMemory<Vector<number> >& mem);
@@ -145,6 +168,9 @@ class BlockMatrixArray : public Subscriptor
                                      * function to call for an object
                                      * created by the default
                                      * constructor.
+                                     *
+                                     * @deprecated the last argument
+                                     * is ignored.
                                      */
     void initialize (const unsigned int n_block_rows,
                     const unsigned int n_block_cols,
@@ -182,12 +208,15 @@ 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().
+                                     *
+                                     * @deprecated the first argument
+                                     * is ignored.
                                      */
     template <class MATRIX>
     void enter_aux (VectorMemory<Vector<number> >& mem,
@@ -308,19 +337,6 @@ class BlockMatrixArray : public Subscriptor
               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
@@ -390,12 +406,6 @@ class BlockMatrixArray : public Subscriptor
                                      * number of blocks per row.
                                      */
     unsigned int block_cols;
-  protected:
-                                    /**
-                                     * The memory used for auxiliary
-                                     * vectors.
-                                     */
-    mutable SmartPointer<VectorMemory<Vector<number> >,BlockMatrixArray<number> > mem;
 };
 
 /*@}*/
@@ -476,6 +486,9 @@ class BlockTrianglePrecondition
                                      * block-quadratic. The additional
                                      * parameter allows for backward
                                      * insertion instead of forward.
+                                     *
+                                     * @deprecated the second argument
+                                     * is ignored.
                                      */
     BlockTrianglePrecondition (unsigned int n_block_rows,
                               VectorMemory<Vector<number> >& mem,
@@ -487,6 +500,9 @@ class BlockTrianglePrecondition
                                      * function to call for an object
                                      * created by the default
                                      * constructor.
+                                     *
+                                     * @deprecated the second argument
+                                     * is ignored.
                                      */
     void initialize (const unsigned int n_block_rows,
                     VectorMemory<Vector<number> >& mem,
@@ -519,6 +535,10 @@ class BlockTrianglePrecondition
                                      * that the diagonal blocks
                                      * should actually be inverse
                                      * matrices or preconditioners.
+                                     *
+                                     * @deprecated The first
+                                     * argument is ignored. User
+                                     * enter() instead.
                                      */
     template <class MATRIX>
     void enter_aux (VectorMemory<Vector<double> >& mem,
@@ -637,26 +657,6 @@ BlockMatrixArray<number>::Entry::Entry (
 
 
 
-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()))
-{}
-
-
-
 template <typename number>
 template <class MATRIX>
 inline
@@ -668,7 +668,6 @@ BlockMatrixArray<number>::enter (
   double prefix,
   bool transpose)
 {
-  Assert (mem != 0, ExcNotInitialized());
   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));
@@ -710,7 +709,6 @@ inline
 void
 BlockMatrixArray<number>::print_latex (STREAM& out) const
 {
-  Assert (mem != 0, ExcNotInitialized());
   out << "\\begin{array}{"
       << std::string(n_block_cols(), 'c')
       << "}" << std::endl;
index 4b779c5052e22177fbe36c97bf767c2128ebb007..6e81c6bfb333d63cf1599213ac085b24bcca3f4c 100644 (file)
@@ -1,8 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
-//    Version: $Name$
 //
-//    Copyright (C) 2005, 2006, 2008 by the deal.II authors
+//    Copyright (C) 2005, 2006, 2008, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -45,20 +44,27 @@ BlockMatrixArray<number>::Entry::~Entry ()
 template <typename number>
 BlockMatrixArray<number>::BlockMatrixArray ()
                : block_rows (0),
-                 block_cols (0),
-                 mem(0, typeid(*this).name())
+                 block_cols (0)
 {}
 
 
 
+template <typename number>
+BlockMatrixArray<number>::BlockMatrixArray (
+  const unsigned int n_block_rows,
+  const unsigned int n_block_cols)
+               : block_rows (n_block_rows),
+                 block_cols (n_block_cols)
+{}
+
+
 template <typename number>
 BlockMatrixArray<number>::BlockMatrixArray (
   const unsigned int n_block_rows,
   const unsigned int n_block_cols,
-  VectorMemory<Vector<number> >& mem)
+  VectorMemory<Vector<number> >&)
                : block_rows (n_block_rows),
-                 block_cols (n_block_cols),
-                 mem(&mem, typeid(*this).name())
+                 block_cols (n_block_cols)
 {}
 
 
@@ -67,14 +73,23 @@ void
 BlockMatrixArray<number>::initialize (
   const unsigned int n_block_rows,
   const unsigned int n_block_cols,
-  VectorMemory<Vector<number> >& memory)
+  VectorMemory<Vector<number> >&)
 {
   block_rows = n_block_rows;
   block_cols = n_block_cols;
-  mem = &memory;
 }
 
 
+template <typename number>
+void
+BlockMatrixArray<number>::initialize (
+  const unsigned int n_block_rows,
+  const unsigned int n_block_cols)
+{
+  block_rows = n_block_rows;
+  block_cols = n_block_cols;
+}
+
 
 
 template <typename number>
@@ -83,7 +98,6 @@ BlockMatrixArray<number>::reinit (
   const unsigned int n_block_rows,
   const unsigned int n_block_cols)
 {
-  Assert (mem != 0, ExcNotInitialized());
   clear();
   block_rows = n_block_rows;
   block_cols = n_block_cols;
@@ -104,14 +118,13 @@ void
 BlockMatrixArray<number>::vmult_add (BlockVector<number>& dst,
                                     const BlockVector<number>& src) const
 {
-  Assert (mem != 0, ExcNotInitialized());
-  
+  GrowingVectorMemory<Vector<number> > mem;
   Assert (dst.n_blocks() == block_rows,
          ExcDimensionMismatch(dst.n_blocks(), block_rows));
   Assert (src.n_blocks() == block_cols,
          ExcDimensionMismatch(src.n_blocks(), block_cols));
 
-  Vector<number>* p_aux = mem->alloc();
+  typename VectorMemory<Vector<number> >::Pointer p_aux(mem);
   Vector<number>& aux = *p_aux;
   
   typename std::vector<Entry>::const_iterator m = entries.begin();
@@ -126,7 +139,6 @@ BlockMatrixArray<number>::vmult_add (BlockVector<number>& dst,
        m->matrix->vmult(aux, src.block(m->col));
       dst.block(m->row).add (m->prefix, aux);
     }
-  mem->free(p_aux);
 }
 
 
@@ -149,7 +161,7 @@ void
 BlockMatrixArray<number>::Tvmult_add (BlockVector<number>& dst,
                                      const BlockVector<number>& src) const
 {
-  Assert (mem != 0, ExcNotInitialized());
+  GrowingVectorMemory<Vector<number> > mem;
   Assert (dst.n_blocks() == block_cols,
          ExcDimensionMismatch(dst.n_blocks(), block_cols));
   Assert (src.n_blocks() == block_rows,
@@ -158,7 +170,7 @@ BlockMatrixArray<number>::Tvmult_add (BlockVector<number>& dst,
   typename std::vector<Entry>::const_iterator m = entries.begin();
   typename std::vector<Entry>::const_iterator end = entries.end();
   
-  Vector<number>* p_aux = mem->alloc();
+  typename VectorMemory<Vector<number> >::Pointer p_aux(mem);
   Vector<number>& aux = *p_aux;
   
   for (; m != end ; ++m)
@@ -170,7 +182,6 @@ BlockMatrixArray<number>::Tvmult_add (BlockVector<number>& dst,
        m->matrix->Tvmult(aux, src.block(m->row));
       dst.block(m->col).add (m->prefix, aux);
     }
-  mem->free(p_aux);
 }
 
 
@@ -193,13 +204,13 @@ BlockMatrixArray<number>::matrix_scalar_product (
   const BlockVector<number>& u,
   const BlockVector<number>& v) const
 {
-  Assert (mem != 0, ExcNotInitialized());
+  GrowingVectorMemory<Vector<number> > mem;
   Assert (u.n_blocks() == block_rows,
          ExcDimensionMismatch(u.n_blocks(), block_rows));
   Assert (v.n_blocks() == block_cols,
          ExcDimensionMismatch(v.n_blocks(), block_cols));
 
-  Vector<number>* p_aux = mem->alloc();
+  typename VectorMemory<Vector<number> >::Pointer p_aux(mem);
   Vector<number>& aux = *p_aux;
   
   typename std::vector<Entry>::const_iterator m;
@@ -221,7 +232,6 @@ BlockMatrixArray<number>::matrix_scalar_product (
        }
       result += u.block(i)*aux;
     }
-  mem->free(p_aux);
 
   return result;
 }
@@ -303,6 +313,7 @@ BlockTrianglePrecondition<number>::do_row (
   BlockVector<number>& dst,
   unsigned int row_num) const
 {
+  GrowingVectorMemory<Vector<number> > mem;
   typename std::vector<typename BlockMatrixArray<number>::Entry>::const_iterator
     m = this->entries.begin();
   typename std::vector<typename BlockMatrixArray<number>::Entry>::const_iterator
@@ -310,7 +321,7 @@ BlockTrianglePrecondition<number>::do_row (
   std::vector<typename std::vector<typename BlockMatrixArray<number>::Entry>::const_iterator>
     diagonals;
   
-  Vector<number>* p_aux = this->mem->alloc();
+  typename VectorMemory<Vector<number> >::Pointer p_aux(mem);
   Vector<number>& aux = *p_aux;
   
   aux.reinit(dst.block(row_num), true);
@@ -372,8 +383,6 @@ BlockTrianglePrecondition<number>::do_row (
        }
       dst.block(row_num) = aux;
     }
-  
-  this->mem->free(p_aux);
 }
 
 

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.