From 0eed0fe694e2940756f8bd95ed689ba613713734 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 17 Jul 2017 14:51:23 -0600 Subject: [PATCH] Avoid raw pointers. Use std::unique_ptr instead. This patch changes SparsityPattern and classes that use it. --- .../lac/chunk_sparse_matrix.templates.h | 12 ++-- .../lac/sparse_decomposition.templates.h | 4 +- include/deal.II/lac/sparse_ilu.templates.h | 12 ++-- include/deal.II/lac/sparse_matrix.templates.h | 24 +++++--- include/deal.II/lac/sparsity_pattern.h | 34 +++++------ source/lac/sparsity_pattern.cc | 57 +++++-------------- 6 files changed, 57 insertions(+), 86 deletions(-) diff --git a/include/deal.II/lac/chunk_sparse_matrix.templates.h b/include/deal.II/lac/chunk_sparse_matrix.templates.h index 0dbb249ac9..62cc5d70ab 100644 --- a/include/deal.II/lac/chunk_sparse_matrix.templates.h +++ b/include/deal.II/lac/chunk_sparse_matrix.templates.h @@ -714,8 +714,8 @@ ChunkSparseMatrix::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::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::matrix_norm_square (const Vector &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::const_iterator v_ptr = v.begin(); for (size_type chunk_row=0; chunk_row::matrix_scalar_product (const Vector &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::const_iterator u_ptr = u.begin(); for (size_type chunk_row=0; chunk_row::residual (Vector &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::iterator dst_ptr = dst.begin(); for (size_type chunk_row=0; chunk_row::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); diff --git a/include/deal.II/lac/sparse_ilu.templates.h b/include/deal.II/lac/sparse_ilu.templates.h index 9f72a0677c..1b3d29e7d0 100644 --- a/include/deal.II/lac/sparse_ilu.templates.h +++ b/include/deal.II/lac/sparse_ilu.templates.h @@ -59,8 +59,8 @@ void SparseILU::initialize (const SparseMatrix &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::val.get(); @@ -147,9 +147,9 @@ void SparseILU::vmult (Vector &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::Tvmult (Vector &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 diff --git a/include/deal.II/lac/sparse_matrix.templates.h b/include/deal.II/lac/sparse_matrix.templates.h index edbe78b3eb..d2a3b13f00 100644 --- a/include/deal.II/lac/sparse_matrix.templates.h +++ b/include/deal.II/lac/sparse_matrix.templates.h @@ -612,7 +612,7 @@ SparseMatrix::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::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::vmult (OutVector &dst, , 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::vmult_add (OutVector &dst, , 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::matrix_norm_square (const Vector &v) const (std::bind (&internal::SparseMatrix::matrix_norm_sqr_on_subrange >, 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::matrix_scalar_product (const Vector &u, (std::bind (&internal::SparseMatrix::matrix_scalar_product_on_subrange >, 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::residual (Vector &dst, (std::bind (&internal::SparseMatrix::residual_sqr_on_subrange ,Vector >, 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)), diff --git a/include/deal.II/lac/sparsity_pattern.h b/include/deal.II/lac/sparsity_pattern.h index 7cd00d7a2c..a7bafba57f 100644 --- a/include/deal.II/lac/sparsity_pattern.h +++ b/include/deal.II/lac/sparsity_pattern.h @@ -31,6 +31,7 @@ #endif #include +#include #include #include #include @@ -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 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 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); } diff --git a/source/lac/sparsity_pattern.cc b/source/lac/sparsity_pattern.cc index 3a51a46b96..6609d8337d 100644 --- a/source/lac/sparsity_pattern.cc +++ b/source/lac/sparsity_pattern.cc @@ -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(), std::placeholders::_1, invalid_entry)); // now allocate the respective memory - size_type *new_colnums = new size_type[nonzero_elements]; + std::unique_ptr 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(&rowstart[0]), -- 2.39.5