]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove Matrix::reinit(). Remove two other deprecated reinit() functions.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 10 May 2004 16:08:11 +0000 (16:08 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 10 May 2004 16:08:11 +0000 (16:08 +0000)
git-svn-id: https://svn.dealii.org/trunk@9224 0785d39b-7218-0410-832d-ea1e28bc413d

14 files changed:
deal.II/doc/news/c-4-0.html
deal.II/examples/step-12/step-12.cc
deal.II/lac/include/lac/block_sparse_matrix.h
deal.II/lac/include/lac/block_sparse_matrix.templates.h
deal.II/lac/include/lac/petsc_matrix_base.h
deal.II/lac/include/lac/petsc_parallel_sparse_matrix.h
deal.II/lac/include/lac/petsc_sparse_matrix.h
deal.II/lac/include/lac/sparse_decomposition.h
deal.II/lac/include/lac/sparse_decomposition.templates.h
deal.II/lac/include/lac/sparse_matrix.h
deal.II/lac/include/lac/sparse_matrix.templates.h
deal.II/lac/include/lac/sparse_mic.h
deal.II/lac/include/lac/sparse_mic.templates.h
deal.II/lac/source/petsc_matrix_base.cc

index 893986159ca0029360c6defe459209d30feb1b0d..cf5fce0dfa39a099537651b910d66c4548c651e7 100644 (file)
@@ -32,7 +32,25 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3 style="color:red">Incompatibilities</h3>
 
 <ol>
-  <li> <p> Changed: All the vector and block vector classes had a member
+  <li> <p> Removed: All the matrix classes have functions <code
+       class="member">reinit</code> that are used to resize the
+       matrix. However, the sparse matrix classes had an equally named
+       function without arguments that basically left the size of the matrix
+       unchanged but set all elements of the matrix to zero. It could also be
+       abused to notify the matrix that the sparsity pattern on which it is
+       built has changed, an inherently risky procedure. The no-argument code
+       <class="member">reinit</code> function has therefore been removed to
+       avoid confusion with the <code class="member">reinit</code> functions
+       that take arguments. Instead, a <code class="member">set_zero</code>
+       function has been introduced that simply sets all elements of the
+       matrix to zero. If you want to notify a sparse matrix that its sparsity
+       pattern has changed, use the <code
+       class="member">reinit(SparsityPattern)</code> function.
+       <br> 
+       (WB 2004/05/10)
+       </p>
+
+  <li> <p> Removed: All the vector and block vector classes had a member
        function <code class="member">clear</code> which simply resets all
        values of the vector to zero. It did not change the size of the vector,
        though. This was confusing, since the standard C++ container classes
@@ -50,6 +68,14 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
        (WB 2004/05/10)
        </p>
 
+  <li> <p> Removed: The <code
+       class="member">SparseLUDecomposition::reinit</code> and <code
+       class="member">SparseMIC::reinit</code> functions without
+       argument had been deprecated before, and have now been removed.
+       <br> 
+       (WB 2004/05/10)
+       </p>
+
   <li> <p> Changed: the template parameter of <code
        class="class">MGTransferPrebuilt</code> is now the complete vector
        type, not just the value type of the vector. This allows to operate
index 8ef66c88b7184755c170dc8aa43bfa8bd22e155c..40b388d51eb7c400c89e0786b911b329435cd93b 100644 (file)
@@ -1624,7 +1624,7 @@ void DGMethod<dim>::run ()
                                       // reinit the system matrix, the
                                       // right hand side vector and
                                       // the Timer object.
-      system_matrix.reinit();
+      system_matrix.set_zero();
       right_hand_side = 0;
       assemble_timer.reset();
 
index 5e86624ae03b159407fa7e3863a6f25540fec51a..3f477e99b2cab9d603a9161654dd094ea93cd688 100644 (file)
@@ -339,24 +339,11 @@ class BlockSparseMatrix : public Subscriptor
     operator = (const BlockSparseMatrix<number> &);
 
                                     /**
-                                     * Reinitialize the object but
-                                     * keep to the sparsity pattern
-                                     * previously used.  This may be
-                                     * necessary if you reinitialized
-                                     * the sparsity structure and
-                                     * want to update the size of the
-                                     * matrix. It only calls
-                                     * SparseMatrix::reinit() on the
-                                     * sub-matrices. The size of this
-                                     * matrix is unchanged.
-                                     *
-                                     * If the sparsity pattern has
-                                     * not changed, then the effect
-                                     * of this function is simply to
-                                     * reset all matrix entries to
-                                     * zero.
+                                     * Set all elements of the matrix to
+                                     * zero, but keep the sparsity pattern
+                                     * previously used.
                                      */
-    virtual void reinit ();
+    void set_zero ();
 
                                     /**
                                      * Reinitialize the sparse matrix
index 144445c123490c291c6cf462e833e7047c2e6962..73b80067ac2cb47a692a129b1a0e0de34cf71a87 100644 (file)
@@ -80,11 +80,11 @@ operator = (const BlockSparseMatrix<number> &m)
 
 template <typename number>
 void
-BlockSparseMatrix<number>::reinit ()
+BlockSparseMatrix<number>::set_zero ()
 {
   for (unsigned int r=0; r<rows; ++r)
     for (unsigned int c=0; c<columns; ++c)
-      block(r,c).reinit ();
+      block(r,c).set_zero ();
 }
 
 
index 32cbcb0cf9027d36b50901dee46c838d85b82a55..51c6ae8cfedc7ba8dcfba5d1e91f1c78aee7ab7f 100644 (file)
@@ -319,7 +319,7 @@ namespace PETScWrappers
                                         * retain the sparsity pattern and all
                                         * the other properties of the matrix.
                                         */
-      void reinit ();
+      void set_zero ();
 
                                        /**
                                         * Set the element (<i>i,j</i>)
index 297d16a8578e3cfbcecc2f9d4b17e3d0dfc45f2c..e7335ce7f77e7f78a4c6062e46623f86b9498445 100644 (file)
@@ -129,7 +129,7 @@ namespace PETScWrappers
                                           * function simply calls the respective
                                           * function of the base class.
                                           */
-        void reinit ();
+        void set_zero ();
 
                                          /**
                                           * Throw away the present matrix and
index 83c8c8744224cda2a18d2b13b00bd7b13a4c03c1..17a5a21b7cfe538851b867c832554655492348d6 100644 (file)
@@ -109,14 +109,6 @@ namespace PETScWrappers
                     const std::vector<unsigned int> &row_lengths,
                     const bool                       is_symmetric = false);
 
-                                       /**
-                                        * Set all matrix entries to zero, but
-                                        * retain the sparsity pattern. This
-                                        * function simply calls the respective
-                                        * function of the base class.
-                                        */
-      void reinit ();
-
                                        /**
                                         * Throw away the present matrix and
                                         * generate one that has the same
@@ -165,18 +157,6 @@ namespace PETScWrappers
                       const std::vector<unsigned int> &row_lengths,
                       const bool                       is_symmetric = false);
   };
-
-
-/// @if NoDoc
-// -------- template and inline functions ----------
-
-  inline
-  void
-  SparseMatrix::reinit ()
-  {
-    MatrixBase::reinit ();
-  }
-///@endif  
 }
 
 #endif // DEAL_II_USE_PETSC
index 70b9074bcbd6d8aadc701cc9e972e95d0a3d1487..6dd9dc9e270cf61691101a26ccb64ce7f3f633fe 100644 (file)
@@ -283,14 +283,6 @@ class SparseLUDecomposition : protected SparseMatrix<number>,
     void initialize (const SparseMatrix<somenumber> &matrix,
                     const AdditionalData parameters);
 
-                                    /**
-                                     * This method is deprecated,
-                                     * and left for backward
-                                     * compability. It will be removed
-                                     * in later versions.
-                                     */
-    void reinit ();
-
                                     /**
                                      * This method is deprecated,
                                      * and left for backward
index 9989f977f31152029fadd6a303c5b344e04daddb..589c4b51ac84ae12815c203456c5319bb6a26b16 100644 (file)
@@ -1,5 +1,5 @@
 //---------------------  sparse_decomposition.templates.h  ----------------
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004
 //    by the deal.II authors and Stephen "Cheffo" Kolaroff
 //
 //    This file is subject to QPL and may not be  distributed
@@ -133,21 +133,6 @@ decompose (const SparseMatrix<somenumber> &matrix,
 
 
 
-template <typename number>
-void
-SparseLUDecomposition<number>::reinit ()
-{
-  decomposed = false;
-  if (true)
-    {
-      std::vector<const unsigned int*> tmp;
-      tmp.swap (prebuilt_lower_bound);
-    };
-  SparseMatrix<number>::reinit ();
-}
-
-
-
 template <typename number>
 void SparseLUDecomposition<number>::reinit (const SparsityPattern& sparsity)
 {
index 8c971a37b14d11a6084d67e94cb27b0ca7583649..4c0f979b5af744963669c004c5e58a126216762a 100644 (file)
@@ -492,31 +492,10 @@ class SparseMatrix : public virtual Subscriptor
     SparseMatrix<number>& operator = (const SparseMatrix<number> &);
 
                                     /**
-                                     * Reinitialize the object but
-                                     * keep to the sparsity pattern
-                                     * previously used.  This may be
-                                     * necessary if the sparsity
-                                     * structure has changed and you
-                                     * want to update the size of the
-                                     * matrix.
-                                     *
-                                     * Note that memory is only
-                                     * reallocated if the new size
-                                     * exceeds the old size. If that
-                                     * is not the case, the allocated
-                                     * memory is not reduced. However,
-                                     * if the sparsity structure is
-                                     * empty (i.e. the dimensions are
-                                     * zero), then all memory is
-                                     * freed.
-                                     *
-                                     * If the sparsity pattern has
-                                     * not changed, then the effect
-                                     * of this function is simply to
-                                     * reset all matrix entries to
-                                     * zero.
+                                     * Set all elements to zero, but keep to
+                                     * the sparsity pattern of the matrix.
                                      */
-    virtual void reinit ();
+    virtual void set_zero ();
 
                                     /**
                                      * Reinitialize the sparse matrix
index e24791371eb7c90a8a01e0a7b4595d29cacd6826..16203a8cf986339f267a0d80416e8e0c6908fd0e 100644 (file)
@@ -80,11 +80,11 @@ SparseMatrix<number>::operator = (const SparseMatrix<number> &m)
 
 template <typename number>
 SparseMatrix<number>::SparseMatrix (const SparsityPattern &c) :
-               cols(&c),
+               cols(0),
                val(0),
                max_len(0)
 {
-  reinit();
+  reinit (c);
 }
 
 
@@ -102,43 +102,47 @@ SparseMatrix<number>::~SparseMatrix ()
 
 template <typename number>
 void
-SparseMatrix<number>::reinit ()
+SparseMatrix<number>::set_zero ()
 {
   Assert (cols != 0, ExcMatrixNotInitialized());
   Assert (cols->compressed || cols->empty(), ExcNotCompressed());
 
-  if (cols->empty()) 
+  if (val)
+    std::fill_n (&val[0], cols->n_nonzero_elements(), 0);
+}
+
+
+
+template <typename number>
+void
+SparseMatrix<number>::reinit (const SparsityPattern &sparsity)
+{
+  cols = &sparsity;
+
+  if (cols->empty())
     {
-      if (val) delete[] val;
+      if (val != 0)
+        delete[] val;
       val = 0;
       max_len = 0;
       return;
-    };
+    }
 
   const unsigned int N = cols->n_nonzero_elements();
   if (N > max_len || max_len == 0)
     {
-      if (val) delete[] val;
+      if (val != 0)
+        delete[] val;
       val = new number[N];
       max_len = N;
-    };
+    }
 
-  if (val)
+  if (val != 0)
     std::fill_n (&val[0], N, 0);
 }
 
 
 
-template <typename number>
-void
-SparseMatrix<number>::reinit (const SparsityPattern &sparsity)
-{
-  cols = &sparsity;
-  reinit ();
-}
-
-
-
 template <typename number>
 void
 SparseMatrix<number>::clear ()
@@ -246,7 +250,7 @@ void
 SparseMatrix<number>::copy_from (const FullMatrix<somenumber> &matrix)
 {
                                   // first delete previous content
-  reinit ();
+  set_zero ();
 
                                   // then copy old matrix
   for (unsigned int row=0; row<matrix.m(); ++row)
index ff808520a916ecaa7b2ea78496a1fff768a5b2b2..bbe5ffc26ede4510717b4079620a859c500acd39 100644 (file)
@@ -1,5 +1,5 @@
 //----------------------------  sparse_mic.h  ---------------------------
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004
 //    by the deal.II authors and Stephen "Cheffo" Kolaroff
 //
 //    This file is subject to QPL and may not be  distributed
@@ -78,14 +78,6 @@ class SparseMIC : public SparseLUDecomposition<number>
     typename SparseLUDecomposition<number>::AdditionalData
     AdditionalData;
 
-                                    /**
-                                     * This method is deprecated, and
-                                     * left for backward
-                                     * compability. It will be
-                                     * removed in later versions.
-                                     */
-    void reinit ();
-
                                     /**
                                      * This method is deprecated, and
                                      * left for backward
index a370fa91ba1694ce6bf33256ed9e8eaadeea406f..f51f45300e01775139da0acab275227ad5d6488c 100644 (file)
@@ -1,5 +1,5 @@
 //----------------------------  sparse_mic.templates.h  ---------------------------
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004
 //    by the deal.II authors and Stephen "Cheffo" Kolaroff
 //
 //    This file is subject to QPL and may not be  distributed
@@ -80,31 +80,6 @@ void SparseMIC<number>::initialize (const SparseMatrix<somenumber> &matrix,
 
 
 
-template <typename number>
-void
-SparseMIC<number>::reinit ()
-{
-  if (true)
-    {
-      std::vector<number> tmp;
-      tmp.swap (diag);
-    };
-  if (true)
-    {
-      std::vector<number> tmp;
-      tmp.swap (inv_diag);
-    };
-  if (true)
-    {
-      std::vector<number> tmp;
-      tmp.swap (inner_sums);
-    };
-
-  SparseLUDecomposition<number>::reinit ();
-}
-
-
-
 template <typename number>
 void SparseMIC<number>::reinit (const SparsityPattern& sparsity)
 {
index 1f0f0fdffd3e7377852d05293fbb1b65cae76277..82cf58217f63b0feca2126dbb14848988859cc2e 100644 (file)
@@ -133,7 +133,7 @@ namespace PETScWrappers
 
 
   void
-  MatrixBase::reinit ()
+  MatrixBase::set_zero ()
   {
     const int ierr = MatZeroEntries (matrix);
     AssertThrow (ierr == 0, ExcPETScError(ierr));    

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.