From 81599415327368d6252b8c683c1337d621724b7d Mon Sep 17 00:00:00 2001 From: David Wells Date: Tue, 30 Jan 2024 17:32:44 -0500 Subject: [PATCH] Use AlignedVector in ReadWriteVector. Why not? --- include/deal.II/lac/read_write_vector.h | 27 ++----- .../deal.II/lac/read_write_vector.templates.h | 77 ++++--------------- 2 files changed, 23 insertions(+), 81 deletions(-) diff --git a/include/deal.II/lac/read_write_vector.h b/include/deal.II/lac/read_write_vector.h index 8effad5f4e..13603f9918 100644 --- a/include/deal.II/lac/read_write_vector.h +++ b/include/deal.II/lac/read_write_vector.h @@ -18,6 +18,7 @@ #include +#include #include #include #include @@ -775,12 +776,6 @@ namespace LinearAlgebra unsigned int global_to_local(const types::global_dof_index global_index) const; - /** - * A helper function that is used to resize the val array. - */ - void - resize_val(const size_type new_allocated_size); - #ifdef DEAL_II_WITH_TRILINOS # ifdef DEAL_II_TRILINOS_WITH_TPETRA /** @@ -818,10 +813,9 @@ namespace LinearAlgebra std::shared_ptr comm_pattern; /** - * Pointer to the array of local elements of this vector. + * Locally stored elements. */ - std::unique_ptr - values; + AlignedVector values; /** * For parallel loops with TBB, this member variable stores the affinity @@ -878,8 +872,6 @@ namespace LinearAlgebra template inline ReadWriteVector::ReadWriteVector() - : Subscriptor() - , values(nullptr, free) { // virtual functions called in constructors and destructors never use the // override in a derived class @@ -893,7 +885,6 @@ namespace LinearAlgebra inline ReadWriteVector::ReadWriteVector( const ReadWriteVector &v) : Subscriptor() - , values(nullptr, free) { this->operator=(v); } @@ -902,8 +893,6 @@ namespace LinearAlgebra template inline ReadWriteVector::ReadWriteVector(const size_type size) - : Subscriptor() - , values(nullptr, free) { // virtual functions called in constructors and destructors never use the // override in a derived class @@ -916,8 +905,6 @@ namespace LinearAlgebra template inline ReadWriteVector::ReadWriteVector( const IndexSet &locally_stored_indices) - : Subscriptor() - , values(nullptr, free) { // virtual functions called in constructors and destructors never use the // override in a derived class @@ -958,7 +945,7 @@ namespace LinearAlgebra inline typename ReadWriteVector::iterator ReadWriteVector::begin() { - return values.get(); + return values.begin(); } @@ -967,7 +954,7 @@ namespace LinearAlgebra inline typename ReadWriteVector::const_iterator ReadWriteVector::begin() const { - return values.get(); + return values.begin(); } @@ -976,7 +963,7 @@ namespace LinearAlgebra inline typename ReadWriteVector::iterator ReadWriteVector::end() { - return values.get() + this->locally_owned_size(); + return values.end(); } @@ -985,7 +972,7 @@ namespace LinearAlgebra inline typename ReadWriteVector::const_iterator ReadWriteVector::end() const { - return values.get() + this->locally_owned_size(); + return values.end(); } diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h index 652990854b..ccc0d9637b 100644 --- a/include/deal.II/lac/read_write_vector.templates.h +++ b/include/deal.II/lac/read_write_vector.templates.h @@ -210,52 +210,13 @@ namespace LinearAlgebra } // namespace internal - template - void - ReadWriteVector::resize_val(const size_type new_alloc_size) - { - if (new_alloc_size == 0) - { - values.reset(); - thread_loop_partitioner = - std::make_shared(); - } - else - { - Number *new_values; - Utilities::System::posix_memalign(reinterpret_cast( - &new_values), - 64, - sizeof(Number) * new_alloc_size); - values.reset(new_values); - - if (new_alloc_size >= 4 * dealii::internal::VectorImplementation:: - minimum_parallel_grain_size) - thread_loop_partitioner = - std::make_shared(); - } - } - - template void ReadWriteVector::reinit(const size_type size, const bool omit_zeroing_entries) { - // check whether we need to reallocate - resize_val(size); - - stored_elements = complete_index_set(size); - stored_elements.compress(); - - // set entries to zero if so requested - if (omit_zeroing_entries == false) - this->operator=(Number()); - - // reset the communication pattern - source_stored_elements.clear(); - comm_pattern.reset(); + reinit(complete_index_set(size), omit_zeroing_entries); } @@ -266,7 +227,8 @@ namespace LinearAlgebra ReadWriteVector::reinit(const ReadWriteVector &v, const bool omit_zeroing_entries) { - resize_val(v.locally_owned_size()); + thread_loop_partitioner = v.thread_loop_partitioner; + values.resize(v.locally_owned_size()); stored_elements = v.get_stored_elements(); @@ -285,10 +247,10 @@ namespace LinearAlgebra ReadWriteVector::reinit(const IndexSet &locally_stored_indices, const bool omit_zeroing_entries) { + thread_loop_partitioner = + std::make_shared(); stored_elements = locally_stored_indices; - - // set vector size and allocate memory - resize_val(stored_elements.n_elements()); + values.resize(stored_elements.n_elements()); // initialize to zero if (omit_zeroing_entries == false) @@ -311,10 +273,7 @@ namespace LinearAlgebra // trilinos data but only if Number=double. Also update documentation that // the argument's lifetime needs to be longer then. If we do this, we need // to think about whether the view should be read/write. - - stored_elements = IndexSet(trilinos_vec.trilinos_partitioner()); - - resize_val(stored_elements.n_elements()); + reinit(IndexSet(trilinos_vec.trilinos_partitioner()), true); TrilinosScalar *start_ptr; int leading_dimension; @@ -322,11 +281,7 @@ namespace LinearAlgebra &leading_dimension); AssertThrow(ierr == 0, ExcTrilinosError(ierr)); - std::copy(start_ptr, start_ptr + leading_dimension, values.get()); - - // reset the communication pattern - source_stored_elements.clear(); - comm_pattern.reset(); + std::copy(start_ptr, start_ptr + leading_dimension, values.data()); } #endif @@ -360,7 +315,7 @@ namespace LinearAlgebra if (locally_owned_size() > 0) { dealii::internal::VectorOperations::Vector_copy copier( - in_vector.values.get(), values.get()); + in_vector.values.data(), values.data()); dealii::internal::VectorOperations::parallel_for( copier, 0, locally_owned_size(), thread_loop_partitioner); } @@ -382,7 +337,7 @@ namespace LinearAlgebra if (locally_owned_size() > 0) { dealii::internal::VectorOperations::Vector_copy copier( - in_vector.values.get(), values.get()); + in_vector.values.data(), values.data()); dealii::internal::VectorOperations::parallel_for( copier, 0, locally_owned_size(), thread_loop_partitioner); } @@ -404,7 +359,7 @@ namespace LinearAlgebra if (this_size > 0) { dealii::internal::VectorOperations::Vector_set setter( - Number(), values.get()); + Number(), values.data()); dealii::internal::VectorOperations::parallel_for( setter, 0, this_size, thread_loop_partitioner); } @@ -739,7 +694,7 @@ namespace LinearAlgebra const double *new_values = target_vector.Values(); const int size = target_vector.MyLength(); - Assert(size == 0 || values != nullptr, + Assert(size == 0 || values.size() > 0, ExcInternalError("Import failed.")); for (int i = 0; i < size; ++i) @@ -754,7 +709,7 @@ namespace LinearAlgebra const double *new_values = target_vector.Values(); const int size = target_vector.MyLength(); - Assert(size == 0 || values != nullptr, + Assert(size == 0 || values.size() > 0, ExcInternalError("Import failed.")); for (int i = 0; i < size; ++i) @@ -769,7 +724,7 @@ namespace LinearAlgebra const double *new_values = target_vector.Values(); const int size = target_vector.MyLength(); - Assert(size == 0 || values != nullptr, + Assert(size == 0 || values.size() > 0, ExcInternalError("Import failed.")); // To ensure that this code also compiles with complex @@ -799,7 +754,7 @@ namespace LinearAlgebra const double *new_values = target_vector.Values(); const int size = target_vector.MyLength(); - Assert(size == 0 || values != nullptr, + Assert(size == 0 || values.size() > 0, ExcInternalError("Import failed.")); for (int i = 0; i < size; ++i) @@ -895,7 +850,7 @@ namespace LinearAlgebra const unsigned int n_elements = stored_elements.n_elements(); if (operation == VectorOperation::insert) { - cudaError_t error_code = cudaMemcpy(values.get(), + cudaError_t error_code = cudaMemcpy(values.data(), cuda_vec.get_values(), n_elements * sizeof(Number), cudaMemcpyDeviceToHost); -- 2.39.5