]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Avoid raw pointers. Use std::unique_ptr instead. 4620/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 17 Jul 2017 20:51:23 +0000 (14:51 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Mon, 17 Jul 2017 20:51:23 +0000 (14:51 -0600)
This patch changes SparsityPattern and classes that use it.

include/deal.II/lac/chunk_sparse_matrix.templates.h
include/deal.II/lac/sparse_decomposition.templates.h
include/deal.II/lac/sparse_ilu.templates.h
include/deal.II/lac/sparse_matrix.templates.h
include/deal.II/lac/sparsity_pattern.h
source/lac/sparsity_pattern.cc

index 0dbb249ac91b7705edd21cac9a98015df0aa666b..62cc5d70aba122220cb912824b88116c5079170f 100644 (file)
@@ -714,8 +714,8 @@ ChunkSparseMatrix<number>::vmult_add (OutVector &dst,
                                            std::cref(*cols),
                                            std::placeholders::_1, std::placeholders::_2,
                                            val.get(),
-                                           cols->sparsity_pattern.rowstart,
-                                           cols->sparsity_pattern.colnums,
+                                           cols->sparsity_pattern.rowstart.get(),
+                                           cols->sparsity_pattern.colnums.get(),
                                            std::cref(src),
                                            std::ref(dst)),
                                 internal::SparseMatrix::minimum_parallel_grain_size/cols->chunk_size+1);
@@ -751,7 +751,7 @@ ChunkSparseMatrix<number>::Tvmult_add (OutVector &dst,
   // like in vmult_add, but don't keep an iterator into dst around since we're
   // not traversing it sequentially this time
   const number    *val_ptr    = val.get();
-  const size_type *colnum_ptr = cols->sparsity_pattern.colnums;
+  const size_type *colnum_ptr = cols->sparsity_pattern.colnums.get();
 
   for (size_type chunk_row=0; chunk_row<n_regular_chunk_rows; ++chunk_row)
     {
@@ -845,7 +845,7 @@ ChunkSparseMatrix<number>::matrix_norm_square (const Vector<somenumber> &v) cons
        n_chunk_rows);
 
   const number    *val_ptr    = val.get();
-  const size_type *colnum_ptr = cols->sparsity_pattern.colnums;
+  const size_type *colnum_ptr = cols->sparsity_pattern.colnums.get();
   typename Vector<somenumber>::const_iterator v_ptr = v.begin();
 
   for (size_type chunk_row=0; chunk_row<n_regular_chunk_rows; ++chunk_row)
@@ -953,7 +953,7 @@ ChunkSparseMatrix<number>::matrix_scalar_product (const Vector<somenumber> &u,
        n_chunk_rows);
 
   const number    *val_ptr    = val.get();
-  const size_type *colnum_ptr = cols->sparsity_pattern.colnums;
+  const size_type *colnum_ptr = cols->sparsity_pattern.colnums.get();
   typename Vector<somenumber>::const_iterator u_ptr = u.begin();
 
   for (size_type chunk_row=0; chunk_row<n_regular_chunk_rows; ++chunk_row)
@@ -1155,7 +1155,7 @@ ChunkSparseMatrix<number>::residual (Vector<somenumber>       &dst,
        n_chunk_rows);
 
   const number    *val_ptr    = val.get();
-  const size_type *colnum_ptr = cols->sparsity_pattern.colnums;
+  const size_type *colnum_ptr = cols->sparsity_pattern.colnums.get();
   typename Vector<somenumber>::iterator dst_ptr = dst.begin();
 
   for (size_type chunk_row=0; chunk_row<n_regular_chunk_rows; ++chunk_row)
index 421d1e0d88853e30368419d90e72ea59a81cd22f..1adc70d93c6d6b5ff56b46950e58ada00712db26 100644 (file)
@@ -135,9 +135,9 @@ void
 SparseLUDecomposition<number>::prebuild_lower_bound()
 {
   const size_type *const
-  column_numbers = this->get_sparsity_pattern().colnums;
+  column_numbers = this->get_sparsity_pattern().colnums.get();
   const std::size_t *const
-  rowstart_indices = this->get_sparsity_pattern().rowstart;
+  rowstart_indices = this->get_sparsity_pattern().rowstart.get();
   const size_type N = this->m();
 
   prebuilt_lower_bound.resize (N);
index 9f72a0677caf35cf50b1e8594f82720b2333fdc6..1b3d29e7d092e20e6a5e9b3de803cb2fa9a1c340 100644 (file)
@@ -59,8 +59,8 @@ void SparseILU<number>::initialize (const SparseMatrix<somenumber> &matrix,
   // translating in essence the algorithm given at the end of section 10.3.2,
   // using the names of variables used there
   const SparsityPattern     &sparsity = this->get_sparsity_pattern();
-  const std::size_t *const ia    = sparsity.rowstart;
-  const size_type *const ja      = sparsity.colnums;
+  const std::size_t *const ia    = sparsity.rowstart.get();
+  const size_type *const ja      = sparsity.colnums.get();
 
   number *luval = this->SparseMatrix<number>::val.get();
 
@@ -147,9 +147,9 @@ void SparseILU<number>::vmult (Vector<somenumber>       &dst,
 
   const size_type N=dst.size();
   const std::size_t *const rowstart_indices
-    = this->get_sparsity_pattern().rowstart;
+    = this->get_sparsity_pattern().rowstart.get();
   const size_type *const column_numbers
-    = this->get_sparsity_pattern().colnums;
+    = this->get_sparsity_pattern().colnums.get();
 
   // solve LUx=b in two steps:
   // first Ly = b, then
@@ -221,9 +221,9 @@ void SparseILU<number>::Tvmult (Vector<somenumber>       &dst,
 
   const size_type N=dst.size();
   const std::size_t *const rowstart_indices
-    = this->get_sparsity_pattern().rowstart;
+    = this->get_sparsity_pattern().rowstart.get();
   const size_type *const column_numbers
-    = this->get_sparsity_pattern().colnums;
+    = this->get_sparsity_pattern().colnums.get();
 
   // solve (LU)'x=b in two steps:
   // first U'y = b, then
index edbe78b3eb7dfb8181a711d57d8fa9e232d41e50..d2a3b13f00808d092299516e3f2e4d19a1f7190c 100644 (file)
@@ -612,7 +612,7 @@ SparseMatrix<number>::add (const size_type  row,
   // unsorted case: first, search all the
   // indices to find out which values we
   // actually need to add.
-  const size_type *const my_cols = cols->colnums;
+  const size_type *const my_cols = cols->colnums.get();
   size_type index = cols->rowstart[row];
   const size_type next_row_index = cols->rowstart[row+1];
 
@@ -671,7 +671,7 @@ SparseMatrix<number>::set (const size_type  row,
   // First, search all the indices to find
   // out which values we actually need to
   // set.
-  const size_type *my_cols = cols->colnums;
+  const size_type *my_cols = cols->colnums.get();
   std::size_t index = cols->rowstart[row], next_index = index;
   const std::size_t next_row_index = cols->rowstart[row+1];
 
@@ -757,8 +757,8 @@ SparseMatrix<number>::vmult (OutVector &dst,
                                            <number,InVector,OutVector>,
                                            std::placeholders::_1, std::placeholders::_2,
                                            val.get(),
-                                           cols->rowstart,
-                                           cols->colnums,
+                                           cols->rowstart.get(),
+                                           cols->colnums.get(),
                                            std::cref(src),
                                            std::ref(dst),
                                            false),
@@ -812,8 +812,8 @@ SparseMatrix<number>::vmult_add (OutVector &dst,
                                            <number,InVector,OutVector>,
                                            std::placeholders::_1, std::placeholders::_2,
                                            val.get(),
-                                           cols->rowstart,
-                                           cols->colnums,
+                                           cols->rowstart.get(),
+                                           cols->colnums.get(),
                                            std::cref(src),
                                            std::ref(dst),
                                            true),
@@ -897,7 +897,9 @@ SparseMatrix<number>::matrix_norm_square (const Vector<somenumber> &v) const
     (std::bind (&internal::SparseMatrix::matrix_norm_sqr_on_subrange
                 <number,Vector<somenumber> >,
                 std::placeholders::_1, std::placeholders::_2,
-                val.get(), cols->rowstart, cols->colnums,
+                val.get(),
+                cols->rowstart.get(),
+                cols->colnums.get(),
                 std::cref(v)),
      0, m(),
      internal::SparseMatrix::minimum_parallel_grain_size);
@@ -960,7 +962,9 @@ SparseMatrix<number>::matrix_scalar_product (const Vector<somenumber> &u,
     (std::bind (&internal::SparseMatrix::matrix_scalar_product_on_subrange
                 <number,Vector<somenumber> >,
                 std::placeholders::_1, std::placeholders::_2,
-                val.get(), cols->rowstart, cols->colnums,
+                val.get(),
+                cols->rowstart.get(),
+                cols->colnums.get(),
                 std::cref(u),
                 std::cref(v)),
      0, m(),
@@ -1347,7 +1351,9 @@ SparseMatrix<number>::residual (Vector<somenumber>       &dst,
                (std::bind (&internal::SparseMatrix::residual_sqr_on_subrange
                            <number,Vector<somenumber>,Vector<somenumber> >,
                            std::placeholders::_1, std::placeholders::_2,
-                           val.get(), cols->rowstart, cols->colnums,
+                           val.get(),
+                           cols->rowstart.get(),
+                           cols->colnums.get(),
                            std::cref(u),
                            std::cref(b),
                            std::ref(dst)),
index 7cd00d7a2c5adbc7002ba7d47bdca368cfcce767..a7bafba57fb7d43fafd68c56bb559db075e1320a 100644 (file)
@@ -31,6 +31,7 @@
 #endif
 #include <boost/serialization/split_member.hpp>
 
+#include <memory>
 #include <vector>
 #include <iostream>
 #include <algorithm>
@@ -479,13 +480,13 @@ public:
    * compressed after this function finishes.
    */
   SparsityPattern (const SparsityPattern  &original,
-                   const unsigned int        max_per_row,
-                   const size_type        extra_off_diagonals);
+                   const unsigned int      max_per_row,
+                   const size_type         extra_off_diagonals);
 
   /**
    * Destructor.
    */
-  ~SparsityPattern ();
+  ~SparsityPattern () = default;
 
   /**
    * Copy operator. For this the same holds as for the copy constructor: it is
@@ -1064,7 +1065,7 @@ private:
    * region that is used. The actual number of elements that was allocated is
    * stored in #max_dim.
    */
-  std::size_t *rowstart;
+  std::unique_ptr<std::size_t[]> rowstart;
 
   /**
    * Array of column numbers. In this array, we store for each non-zero
@@ -1088,7 +1089,7 @@ private:
    * sorted, such that finding whether an element exists and determining its
    * position can be done by a binary search.
    */
-  size_type *colnums;
+  std::unique_ptr<size_type[]> colnums;
 
   /**
    * Store whether the compress() function was called for this object.
@@ -1164,10 +1165,10 @@ namespace SparsityPatternIterators
     Assert (is_valid_entry() == true, ExcInvalidIterator());
 
     const std::size_t *insert_point =
-      std::upper_bound(sparsity_pattern->rowstart,
-                       sparsity_pattern->rowstart + sparsity_pattern->rows + 1,
+      std::upper_bound(sparsity_pattern->rowstart.get(),
+                       sparsity_pattern->rowstart.get() + sparsity_pattern->rows + 1,
                        index_within_sparsity);
-    return insert_point - sparsity_pattern->rowstart - 1;
+    return insert_point - sparsity_pattern->rowstart.get() - 1;
   }
 
 
@@ -1425,8 +1426,8 @@ SparsityPattern::save (Archive &ar, const unsigned int) const
 
   ar &max_dim &rows &cols &max_vec_len &max_row_length &compressed &store_diagonal_first_in_row;
 
-  ar &boost::serialization::make_array(rowstart, max_dim + 1);
-  ar &boost::serialization::make_array(colnums, max_vec_len);
+  ar &boost::serialization::make_array(rowstart.get(), max_dim + 1);
+  ar &boost::serialization::make_array(colnums.get(), max_vec_len);
 }
 
 
@@ -1441,16 +1442,11 @@ SparsityPattern::load (Archive &ar, const unsigned int)
 
   ar &max_dim &rows &cols &max_vec_len &max_row_length &compressed &store_diagonal_first_in_row;
 
-  if (rowstart != nullptr)
-    delete[] rowstart;
-  rowstart = new std::size_t[max_dim + 1];
+  rowstart.reset (new std::size_t[max_dim + 1]);
+  colnums.reset (new size_type[max_vec_len]);
 
-  if (colnums != nullptr)
-    delete[] colnums;
-  colnums = new size_type[max_vec_len];
-
-  ar &boost::serialization::make_array(rowstart, max_dim + 1);
-  ar &boost::serialization::make_array(colnums, max_vec_len);
+  ar &boost::serialization::make_array(rowstart.get(), max_dim + 1);
+  ar &boost::serialization::make_array(colnums.get(), max_vec_len);
 }
 
 
index 3a51a46b9666f29a3cda20186054f7a8ccec22de..6609d8337d43341b65b4660bf25c06b12cb90f61 100644 (file)
@@ -214,14 +214,6 @@ SparsityPattern::SparsityPattern (const SparsityPattern &original,
 
 
 
-SparsityPattern::~SparsityPattern ()
-{
-  if (rowstart != nullptr)  delete[] rowstart;
-  if (colnums != nullptr)   delete[] colnums;
-}
-
-
-
 SparsityPattern &
 SparsityPattern::operator = (const SparsityPattern &s)
 {
@@ -265,14 +257,14 @@ SparsityPattern::reinit (const size_type m,
   // delete empty matrices
   if ((m==0) || (n==0))
     {
-      if (rowstart)  delete[] rowstart;
-      if (colnums)   delete[] colnums;
-      rowstart = nullptr;
-      colnums = nullptr;
+      rowstart.reset();
+      colnums.reset();
+
       max_vec_len = max_dim = rows = cols = 0;
       // if dimension is zero: ignore max_per_row
       max_row_length = 0;
       compressed = false;
+
       return;
     }
 
@@ -301,14 +293,8 @@ SparsityPattern::reinit (const size_type m,
   if (vec_len == 0)
     {
       vec_len = 1;
-      if (colnums)
-        {
-          delete[] colnums;
-          colnums = nullptr;
-        }
-
       max_vec_len = vec_len;
-      colnums = new size_type[max_vec_len];
+      colnums.reset (new size_type[max_vec_len]);
     }
 
   max_row_length = (row_lengths.size() == 0 ?
@@ -328,27 +314,15 @@ SparsityPattern::reinit (const size_type m,
   // and try to delete the memory a second time.
   if (rows > max_dim)
     {
-      if (rowstart)
-        {
-          delete[] rowstart;
-          rowstart = nullptr;
-        }
-
       max_dim = rows;
-      rowstart = new std::size_t[max_dim+1];
+      rowstart.reset (new std::size_t[max_dim+1]);
     }
 
   // allocate memory for the column numbers if necessary
   if (vec_len > max_vec_len)
     {
-      if (colnums)
-        {
-          delete[] colnums;
-          colnums = nullptr;
-        }
-
       max_vec_len = vec_len;
-      colnums = new size_type[max_vec_len];
+      colnums.reset (new size_type[max_vec_len]);
     }
 
   // set the rowstart array
@@ -398,7 +372,7 @@ SparsityPattern::compress ()
                      &colnums[rowstart[rows]],
                      std::bind(std::not_equal_to<size_type>(), std::placeholders::_1, invalid_entry));
   // now allocate the respective memory
-  size_type *new_colnums = new size_type[nonzero_elements];
+  std::unique_ptr<size_type[]> new_colnums (new size_type[nonzero_elements]);
 
 
   // reserve temporary storage to store the entries of one row
@@ -470,9 +444,9 @@ SparsityPattern::compress ()
   // set iterator-past-the-end
   rowstart[rows] = next_row_start;
 
-  // set colnums to the newly allocated array and delete the old one
-  delete[] colnums;
-  colnums = new_colnums;
+  // set colnums to the newly allocated array and delete previous content
+  // in the process
+  colnums = std::move(new_colnums);
 
   // store the size
   max_vec_len = nonzero_elements;
@@ -1000,13 +974,8 @@ SparsityPattern::block_read (std::istream &in)
   AssertThrow (c == '[', ExcIO());
 
   // reallocate space
-  if (rowstart)
-    delete[] rowstart;
-  if (colnums)
-    delete[] colnums;
-
-  rowstart = new std::size_t[max_dim+1];
-  colnums  = new size_type[max_vec_len];
+  rowstart.reset (new std::size_t[max_dim+1]);
+  colnums.reset (new size_type[max_vec_len]);
 
   // then read data
   in.read (reinterpret_cast<char *>(&rowstart[0]),

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.