From d716983f2a91177d5a2364ec6ba5d8c1c8b25175 Mon Sep 17 00:00:00 2001 From: David Wells Date: Sun, 20 Jan 2019 22:32:32 -0500 Subject: [PATCH] Use AlignedVector to manage Vector's memory. Vector already requires that the underlying buffer be aligned and uses exactly the same function (posix_memalign) as AlignedVector to do it: we can simplify the implementation of Vector a lot by just using an AlignedVector to manage the data instead. --- include/deal.II/lac/vector.h | 102 ++----- include/deal.II/lac/vector.templates.h | 383 +++++++++---------------- tests/serialization/vector.output | 4 +- 3 files changed, 165 insertions(+), 324 deletions(-) diff --git a/include/deal.II/lac/vector.h b/include/deal.II/lac/vector.h index d44d28ee3b..2dadf1da5e 100644 --- a/include/deal.II/lac/vector.h +++ b/include/deal.II/lac/vector.h @@ -19,6 +19,7 @@ #include +#include #include #include #include @@ -29,14 +30,6 @@ #include #include -// boost::serialization::make_array used to be in array.hpp, but was -// moved to a different file in BOOST 1.64 -#include -#if BOOST_VERSION >= 106400 -# include -#else -# include -#endif #include #include @@ -152,6 +145,8 @@ public: * * We would like to make this constructor explicit, but standard containers * insist on using it implicitly. + * + * @dealiiOperationIsMultithreaded */ Vector(const Vector &v); @@ -159,19 +154,13 @@ public: * Move constructor. Creates a new vector by stealing the internal data of * the vector @p v. */ - Vector(Vector &&v) noexcept; + Vector(Vector &&v) noexcept = default; /** * Copy constructor taking a vector of another data type. This will fail if * there is no conversion path from @p OtherNumber to @p Number. Note that * you may lose accuracy when copying to a vector with data elements with * less accuracy. - * - * Older versions of gcc did not honor the @p explicit keyword on template - * constructors. In such cases, it is easy to accidentally write code that - * can be very inefficient, since the compiler starts performing hidden - * conversions. To avoid this, this function is disabled if we have detected - * a broken compiler during configuration. */ template explicit Vector(const Vector &v); @@ -365,7 +354,7 @@ public: * have after being newly default-constructed. */ Vector & - operator=(Vector &&v) noexcept; + operator=(Vector &&v) noexcept = default; /** * Copy the given vector. Resize the present vector if necessary. @@ -993,27 +982,9 @@ public: private: /** - * Dimension. Actual number of components contained in the vector. Get this - * number by calling size(). - */ - size_type vec_size; - - /** - * Amount of memory actually reserved for this vector. This number may be - * greater than @p vec_size if a @p reinit was called with less memory - * requirements than the vector needed last time. At present @p reinit does - * not free memory when the number of needed elements is reduced. - */ - size_type max_vec_size; - - /** - * Pointer to the array of elements of this vector. - * - * Because we allocate these arrays via Utilities::System::posix_memalign, - * we need to use a custom deleter for this object that does not call - * delete[], but instead calls @p free(). + * Array of elements owned by this vector. */ - std::unique_ptr values; + AlignedVector values; /** * For parallel loops with TBB, this member variable stores the affinity @@ -1022,25 +993,11 @@ private: mutable std::shared_ptr thread_loop_partitioner; - /** - * Allocate and align @p values along 64-byte boundaries. The size of the - * allocated memory is determined by @p max_vec_size . Copy first - * @p copy_n_el from the old values. - */ - void - allocate(const size_type copy_n_el = 0); - /** * Make all other vector types friends. */ template friend class Vector; - - /** - * LAPACK matrices need access to the data. - */ - template - friend class LAPACKFullMatrix; }; /*@}*/ @@ -1060,9 +1017,6 @@ Vector::lp_norm(const real_type) const; template inline Vector::Vector() - : vec_size(0) - , max_vec_size(0) - , values(nullptr, &free) { // virtual functions called in constructors and destructors never use the // override in a derived class @@ -1075,9 +1029,6 @@ inline Vector::Vector() template template Vector::Vector(const InputIterator first, const InputIterator last) - : vec_size(0) - , max_vec_size(0) - , values(nullptr, &free) { // allocate memory. do not initialize it, as we will copy over to it in a // second @@ -1089,9 +1040,6 @@ Vector::Vector(const InputIterator first, const InputIterator last) template inline Vector::Vector(const size_type n) - : vec_size(0) - , max_vec_size(0) - , values(nullptr, &free) { // virtual functions called in constructors and destructors never use the // override in a derived class @@ -1105,7 +1053,7 @@ template inline typename Vector::size_type Vector::size() const { - return vec_size; + return values.size(); } @@ -1122,7 +1070,7 @@ template inline typename Vector::pointer Vector::data() { - return values.get(); + return values.data(); } @@ -1131,7 +1079,7 @@ template inline typename Vector::const_pointer Vector::data() const { - return values.get(); + return values.data(); } @@ -1140,7 +1088,7 @@ template inline typename Vector::iterator Vector::begin() { - return values.get(); + return values.begin(); } @@ -1149,7 +1097,7 @@ template inline typename Vector::const_iterator Vector::begin() const { - return values.get(); + return values.begin(); } @@ -1158,7 +1106,7 @@ template inline typename Vector::iterator Vector::end() { - return values.get() + vec_size; + return values.end(); } @@ -1167,7 +1115,7 @@ template inline typename Vector::const_iterator Vector::end() const { - return values.get() + vec_size; + return values.end(); } @@ -1176,7 +1124,7 @@ template inline Number Vector::operator()(const size_type i) const { - Assert(i < vec_size, ExcIndexRange(i, 0, vec_size)); + Assert(i < size(), ExcIndexRange(i, 0, size())); return values[i]; } @@ -1186,7 +1134,7 @@ template inline Number & Vector::operator()(const size_type i) { - Assert(i < vec_size, ExcIndexRangeType(i, 0, vec_size)); + Assert(i < size(), ExcIndexRangeType(i, 0, size())); return values[i]; } @@ -1271,7 +1219,7 @@ Vector::add(const std::vector &indices, { Assert(indices.size() == values.size(), ExcDimensionMismatch(indices.size(), values.size())); - add(indices.size(), indices.data(), values.values.get()); + add(indices.size(), indices.data(), values.values.begin()); } @@ -1285,7 +1233,7 @@ Vector::add(const size_type n_indices, { for (size_type i = 0; i < n_indices; ++i) { - Assert(indices[i] < vec_size, ExcIndexRange(indices[i], 0, vec_size)); + Assert(indices[i] < size(), ExcIndexRange(indices[i], 0, size())); Assert( numbers::is_finite(values[i]), ExcMessage( @@ -1327,8 +1275,6 @@ template inline void Vector::swap(Vector &v) { - std::swap(vec_size, v.vec_size); - std::swap(max_vec_size, v.max_vec_size); std::swap(values, v.values); } @@ -1341,9 +1287,7 @@ Vector::save(Archive &ar, const unsigned int) const { // forward to serialization function in the base class. ar &static_cast(*this); - - ar &vec_size &max_vec_size; - ar & boost::serialization::make_array(values.get(), max_vec_size); + ar &values; } @@ -1353,15 +1297,9 @@ template inline void Vector::load(Archive &ar, const unsigned int) { - // get rid of previous content - values.reset(); - // the load stuff again from the archive ar &static_cast(*this); - ar &vec_size &max_vec_size; - - allocate(); - ar &boost::serialization::make_array(values.get(), max_vec_size); + ar &values; } #endif diff --git a/include/deal.II/lac/vector.templates.h b/include/deal.II/lac/vector.templates.h index e67ffd71ae..87c126c8b8 100644 --- a/include/deal.II/lac/vector.templates.h +++ b/include/deal.II/lac/vector.templates.h @@ -46,15 +46,8 @@ DEAL_II_NAMESPACE_OPEN template Vector::Vector(const Vector &v) : Subscriptor() - , vec_size(v.size()) - , max_vec_size(v.size()) - , values(nullptr, &free) { - if (vec_size != 0) - { - allocate(); - *this = v; - } + *this = v; } @@ -73,34 +66,11 @@ Vector::apply_givens_rotation(const std::array &csr, -template -Vector::Vector(Vector &&v) noexcept - : Subscriptor(std::move(v)) - , vec_size(v.vec_size) - , max_vec_size(v.max_vec_size) - , values(std::move(v.values)) - , thread_loop_partitioner(std::move(v.thread_loop_partitioner)) -{ - v.vec_size = 0; - v.max_vec_size = 0; - v.values = nullptr; -} - - - template template Vector::Vector(const Vector &v) - : Subscriptor() - , vec_size(v.size()) - , max_vec_size(v.size()) - , values(nullptr, &free) { - if (vec_size != 0) - { - allocate(); - *this = v; - } + *this = v; } @@ -154,10 +124,6 @@ namespace internal template Vector::Vector(const PETScWrappers::VectorBase &v) - : Subscriptor() - , vec_size(0) - , max_vec_size(0) - , values(nullptr, &free) { if (v.size() != 0) { @@ -171,15 +137,10 @@ Vector::Vector(const PETScWrappers::VectorBase &v) template Vector::Vector(const TrilinosWrappers::MPI::Vector &v) - : Subscriptor() - , vec_size(v.size()) - , max_vec_size(v.size()) - , values(nullptr, &free) + : values(v.size()) { - if (vec_size != 0) + if (size() != 0) { - allocate(); - // Copy the distributed vector to // a local one at all processors // that know about the original vector. @@ -188,12 +149,12 @@ Vector::Vector(const TrilinosWrappers::MPI::Vector &v) // this, but it has not yet been // found. TrilinosWrappers::MPI::Vector localized_vector; - localized_vector.reinit(complete_index_set(vec_size), + localized_vector.reinit(complete_index_set(size()), v.get_mpi_communicator()); localized_vector.reinit(v, false, true); - Assert(localized_vector.size() == vec_size, - ExcDimensionMismatch(localized_vector.size(), vec_size)); + Assert(localized_vector.size() == size(), + ExcDimensionMismatch(localized_vector.size(), size())); // get a representation of the vector // and copy it @@ -202,7 +163,7 @@ Vector::Vector(const TrilinosWrappers::MPI::Vector &v) int ierr = localized_vector.trilinos_vector().ExtractView(&start_ptr); AssertThrow(ierr == 0, ExcTrilinosError(ierr)); - std::copy(start_ptr[0], start_ptr[0] + vec_size, begin()); + std::copy(start_ptr[0], start_ptr[0] + size(), begin()); } } @@ -217,16 +178,16 @@ Vector::operator=(const Vector &v) return *this; thread_loop_partitioner = v.thread_loop_partitioner; - if (vec_size != v.vec_size) + if (size() != v.size()) reinit(v, true); - if (vec_size > 0) + if (0 < size()) { dealii::internal::VectorOperations::Vector_copy copier( - v.values.get(), values.get()); + v.begin(), begin()); internal::VectorOperations::parallel_for(copier, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -235,39 +196,20 @@ Vector::operator=(const Vector &v) -template -inline Vector & -Vector::operator=(Vector &&v) noexcept -{ - Subscriptor::operator=(std::move(v)); - - vec_size = v.vec_size; - max_vec_size = v.max_vec_size; - values = std::move(v.values); - thread_loop_partitioner = std::move(v.thread_loop_partitioner); - - v.vec_size = 0; - v.max_vec_size = 0; - - return *this; -} - - - template template inline Vector & Vector::operator=(const Vector &v) { thread_loop_partitioner = v.thread_loop_partitioner; - if (vec_size != v.vec_size) + if (size() != v.size()) reinit(v, true); dealii::internal::VectorOperations::Vector_copy copier( - v.values.get(), values.get()); + v.begin(), begin()); internal::VectorOperations::parallel_for(copier, 0, - vec_size, + size(), thread_loop_partitioner); return *this; @@ -279,35 +221,41 @@ template inline void Vector::reinit(const size_type n, const bool omit_zeroing_entries) { - if (n == 0) + const std::size_t old_size = size(); + + // avoid allocating if the new size is not larger than the old size + if (n <= size()) { - values.reset(); - max_vec_size = vec_size = 0; thread_loop_partitioner = std::make_shared(); + if (n == 0) + { + values.clear(); + } + else if (n != size()) + values.resize_fast(n); + + if (!omit_zeroing_entries) + values.fill(); return; } - if (n > max_vec_size) - { - max_vec_size = n; - allocate(); - } + // otherwise size() < n and we must allocate + AlignedVector new_values; + new_values.resize_fast(n); + if (!omit_zeroing_entries) + new_values.fill(); + new_values.swap(values); - if (vec_size != n) + if (old_size != size()) { - vec_size = n; - // only reset the partitioner if we actually expect a significant vector // size - if (vec_size >= + if (size() >= 4 * internal::VectorImplementation::minimum_parallel_grain_size) thread_loop_partitioner = std::make_shared(); } - - if (omit_zeroing_entries == false) - *this = Number(); } @@ -316,37 +264,26 @@ template inline void Vector::grow_or_shrink(const size_type n) { + const std::size_t old_size = size(); if (n == 0) { - values.reset(); - max_vec_size = vec_size = 0; + values.clear(); thread_loop_partitioner = std::make_shared(); return; } - const size_type s = std::min(vec_size, n); - if (n > max_vec_size) - { - max_vec_size = n; - allocate(s); - } + values.resize(n); - if (vec_size != n) + if (old_size != n) { - vec_size = n; - // only reset the partitioner if we actually expect a significant vector // size - if (vec_size >= + if (size() >= 4 * internal::VectorImplementation::minimum_parallel_grain_size) thread_loop_partitioner = std::make_shared(); } - - // pad with zeroes - for (size_type i = s; i < vec_size; ++i) - values[i] = Number(); } @@ -359,21 +296,7 @@ Vector::reinit(const Vector &v, { thread_loop_partitioner = v.thread_loop_partitioner; - if (v.vec_size == 0) - { - values.reset(); - max_vec_size = vec_size = 0; - return; - } - - if (v.vec_size > max_vec_size) - { - max_vec_size = v.vec_size; - allocate(); - } - vec_size = v.vec_size; - if (omit_zeroing_entries == false) - *this = Number(); + reinit(v.size(), omit_zeroing_entries); } @@ -382,9 +305,9 @@ template bool Vector::all_zero() const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) if (values[i] != Number()) return false; return true; @@ -396,9 +319,9 @@ template bool Vector::is_non_negative() const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) if (!internal::VectorOperations::is_non_negative(values[i])) return false; @@ -413,14 +336,14 @@ Vector::operator=(const Number s) { AssertIsFinite(s); if (s != Number()) - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); - if (vec_size > 0) + if (size() > 0) { - internal::VectorOperations::Vector_set setter(s, values.get()); + internal::VectorOperations::Vector_set setter(s, values.begin()); internal::VectorOperations::parallel_for(setter, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -435,14 +358,14 @@ Vector::operator*=(const Number factor) { AssertIsFinite(factor); - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); internal::VectorOperations::Vectorization_multiply_factor - vector_multiply(values.get(), factor); + vector_multiply(values.begin(), factor); internal::VectorOperations::parallel_for(vector_multiply, 0, - vec_size, + size(), thread_loop_partitioner); return *this; @@ -456,14 +379,14 @@ Vector::add(const Number a, const Vector &v) { AssertIsFinite(a); - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == v.vec_size, ExcDimensionMismatch(vec_size, v.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); internal::VectorOperations::Vectorization_add_av vector_add_av( - values.get(), v.values.get(), a); + values.begin(), v.values.begin(), a); internal::VectorOperations::parallel_for(vector_add_av, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -476,14 +399,14 @@ Vector::sadd(const Number x, const Number a, const Vector &v) AssertIsFinite(x); AssertIsFinite(a); - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == v.vec_size, ExcDimensionMismatch(vec_size, v.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); internal::VectorOperations::Vectorization_sadd_xav vector_sadd_xav( - values.get(), v.values.get(), a, x); + values.begin(), v.values.begin(), a, x); internal::VectorOperations::parallel_for(vector_sadd_xav, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -493,18 +416,18 @@ template template Number Vector::operator*(const Vector &v) const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); if (PointerComparison::equal(this, &v)) return norm_sqr(); - Assert(vec_size == v.size(), ExcDimensionMismatch(vec_size, v.size())); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); Number sum; - internal::VectorOperations::Dot dot(values.get(), - v.values.get()); + internal::VectorOperations::Dot dot(values.begin(), + v.values.begin()); internal::VectorOperations::parallel_reduce( - dot, 0, vec_size, sum, thread_loop_partitioner); + dot, 0, size(), sum, thread_loop_partitioner); AssertIsFinite(sum); return sum; @@ -516,12 +439,12 @@ template typename Vector::real_type Vector::norm_sqr() const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); real_type sum; - internal::VectorOperations::Norm2 norm2(values.get()); + internal::VectorOperations::Norm2 norm2(values.begin()); internal::VectorOperations::parallel_reduce( - norm2, 0, vec_size, sum, thread_loop_partitioner); + norm2, 0, size(), sum, thread_loop_partitioner); AssertIsFinite(sum); @@ -534,12 +457,12 @@ template Number Vector::mean_value() const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); Number sum; - internal::VectorOperations::MeanValue mean(values.get()); + internal::VectorOperations::MeanValue mean(values.begin()); internal::VectorOperations::parallel_reduce( - mean, 0, vec_size, sum, thread_loop_partitioner); + mean, 0, size(), sum, thread_loop_partitioner); return sum / real_type(size()); } @@ -550,12 +473,12 @@ template typename Vector::real_type Vector::l1_norm() const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); real_type sum; - internal::VectorOperations::Norm1 norm1(values.get()); + internal::VectorOperations::Norm1 norm1(values.begin()); internal::VectorOperations::parallel_reduce( - norm1, 0, vec_size, sum, thread_loop_partitioner); + norm1, 0, size(), sum, thread_loop_partitioner); return sum; } @@ -571,12 +494,12 @@ Vector::l2_norm() const // might still be finite. In that case, recompute it (this is a rare case, // so working on the vector twice is uncritical and paid off by the extended // precision) using the BLAS approach with a weight, see e.g. dnrm2.f. - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); real_type norm_square; - internal::VectorOperations::Norm2 norm2(values.get()); + internal::VectorOperations::Norm2 norm2(values.begin()); internal::VectorOperations::parallel_reduce( - norm2, 0, vec_size, norm_square, thread_loop_partitioner); + norm2, 0, size(), norm_square, thread_loop_partitioner); if (numbers::is_finite(norm_square) && norm_square >= std::numeric_limits::min()) return static_cast::real_type>( @@ -585,7 +508,7 @@ Vector::l2_norm() const { real_type scale = 0.; real_type sum = 1.; - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) { if (values[i] != Number()) { @@ -612,7 +535,7 @@ template typename Vector::real_type Vector::lp_norm(const real_type p) const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); if (p == 1.) return l1_norm(); @@ -620,9 +543,9 @@ Vector::lp_norm(const real_type p) const return l2_norm(); real_type sum; - internal::VectorOperations::NormP normp(values.get(), p); + internal::VectorOperations::NormP normp(values.begin(), p); internal::VectorOperations::parallel_reduce( - normp, 0, vec_size, sum, thread_loop_partitioner); + normp, 0, size(), sum, thread_loop_partitioner); if (numbers::is_finite(sum) && sum >= std::numeric_limits::min()) return std::pow(sum, static_cast(1. / p)); @@ -630,7 +553,7 @@ Vector::lp_norm(const real_type p) const { real_type scale = 0.; real_type sum = 1.; - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) { if (values[i] != Number()) { @@ -655,11 +578,11 @@ template typename Vector::real_type Vector::linfty_norm() const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); real_type max = 0.; - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) max = std::max(numbers::NumberTraits::abs(values[i]), max); return max; @@ -673,17 +596,17 @@ Vector::add_and_dot(const Number a, const Vector &V, const Vector &W) { - Assert(vec_size != 0, ExcEmptyObject()); - AssertDimension(vec_size, V.size()); - AssertDimension(vec_size, W.size()); + Assert(size() != 0, ExcEmptyObject()); + AssertDimension(size(), V.size()); + AssertDimension(size(), W.size()); Number sum; - internal::VectorOperations::AddAndDot adder(this->values.get(), - V.values.get(), - W.values.get(), + internal::VectorOperations::AddAndDot adder(values.begin(), + V.values.begin(), + W.values.begin(), a); internal::VectorOperations::parallel_reduce( - adder, 0, vec_size, sum, thread_loop_partitioner); + adder, 0, size(), sum, thread_loop_partitioner); AssertIsFinite(sum); return sum; @@ -695,14 +618,14 @@ template Vector & Vector::operator+=(const Vector &v) { - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == v.vec_size, ExcDimensionMismatch(vec_size, v.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); internal::VectorOperations::Vectorization_add_v vector_add( - values.get(), v.values.get()); + values.begin(), v.values.begin()); internal::VectorOperations::parallel_for(vector_add, 0, - vec_size, + size(), thread_loop_partitioner); return *this; } @@ -713,14 +636,14 @@ template Vector & Vector::operator-=(const Vector &v) { - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == v.vec_size, ExcDimensionMismatch(vec_size, v.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); internal::VectorOperations::Vectorization_subtract_v vector_subtract( - values.get(), v.values.get()); + values.begin(), v.values.begin()); internal::VectorOperations::parallel_for(vector_subtract, 0, - vec_size, + size(), thread_loop_partitioner); return *this; @@ -732,13 +655,13 @@ template void Vector::add(const Number v) { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); internal::VectorOperations::Vectorization_add_factor vector_add( - values.get(), v); + values.begin(), v); internal::VectorOperations::parallel_for(vector_add, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -754,15 +677,15 @@ Vector::add(const Number a, AssertIsFinite(a); AssertIsFinite(b); - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == v.vec_size, ExcDimensionMismatch(vec_size, v.vec_size)); - Assert(vec_size == w.vec_size, ExcDimensionMismatch(vec_size, w.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); + Assert(size() == w.size(), ExcDimensionMismatch(size(), w.size())); internal::VectorOperations::Vectorization_add_avpbw vector_add( - values.get(), v.values.get(), w.values.get(), a, b); + values.begin(), v.values.begin(), w.values.begin(), a, b); internal::VectorOperations::parallel_for(vector_add, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -774,14 +697,14 @@ Vector::sadd(const Number x, const Vector &v) { AssertIsFinite(x); - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == v.vec_size, ExcDimensionMismatch(vec_size, v.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); internal::VectorOperations::Vectorization_sadd_xv vector_sadd( - values.get(), v.values.get(), x); + values.begin(), v.values.begin(), x); internal::VectorOperations::parallel_for(vector_sadd, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -791,14 +714,14 @@ template void Vector::scale(const Vector &s) { - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == s.vec_size, ExcDimensionMismatch(vec_size, s.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == s.size(), ExcDimensionMismatch(size(), s.size())); internal::VectorOperations::Vectorization_scale vector_scale( - values.get(), s.values.get()); + values.begin(), s.values.begin()); internal::VectorOperations::parallel_for(vector_scale, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -809,10 +732,10 @@ template void Vector::scale(const Vector &s) { - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == s.vec_size, ExcDimensionMismatch(vec_size, s.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == s.size(), ExcDimensionMismatch(size(), s.size())); - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) values[i] *= Number(s.values[i]); } @@ -824,14 +747,14 @@ Vector::equ(const Number a, const Vector &u) { AssertIsFinite(a); - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == u.vec_size, ExcDimensionMismatch(vec_size, u.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == u.size(), ExcDimensionMismatch(size(), u.size())); internal::VectorOperations::Vectorization_equ_au vector_equ( - values.get(), u.values.get(), a); + values.begin(), u.values.begin(), a); internal::VectorOperations::parallel_for(vector_equ, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -844,8 +767,8 @@ Vector::equ(const Number a, const Vector &u) { AssertIsFinite(a); - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == u.vec_size, ExcDimensionMismatch(vec_size, u.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(size() == u.size(), ExcDimensionMismatch(size(), u.size())); // set the result vector to a*u. we have to // convert the elements of u to the type of @@ -853,7 +776,7 @@ Vector::equ(const Number a, const Vector &u) // because // operator*(complex,complex) // is not defined by default - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) values[i] = a * Number(u.values[i]); } @@ -863,19 +786,18 @@ template void Vector::ratio(const Vector &a, const Vector &b) { - Assert(vec_size != 0, ExcEmptyObject()); - Assert(a.vec_size == b.vec_size, - ExcDimensionMismatch(a.vec_size, b.vec_size)); + Assert(size() != 0, ExcEmptyObject()); + Assert(a.size() == b.size(), ExcDimensionMismatch(a.size(), b.size())); // no need to reinit with zeros, since // we overwrite them anyway reinit(a.size(), true); internal::VectorOperations::Vectorization_ratio vector_ratio( - values.get(), a.values.get(), b.values.get()); + values.begin(), a.begin(), b.begin()); internal::VectorOperations::parallel_for(vector_ratio, 0, - vec_size, + size(), thread_loop_partitioner); } @@ -885,7 +807,7 @@ template Vector & Vector::operator=(const BlockVector &v) { - if (v.size() != vec_size) + if (v.size() != size()) reinit(v.size(), true); size_type this_index = 0; @@ -915,9 +837,9 @@ template Vector & Vector::operator=(const TrilinosWrappers::MPI::Vector &v) { - if (v.size() != vec_size) + if (v.size() != size()) reinit(v.size(), true); - if (vec_size != 0) + if (size() != 0) { // Copy the distributed vector to // a local one at all processors @@ -927,12 +849,12 @@ Vector::operator=(const TrilinosWrappers::MPI::Vector &v) // this, but it has not yet been // found. TrilinosWrappers::MPI::Vector localized_vector; - localized_vector.reinit(complete_index_set(vec_size), + localized_vector.reinit(complete_index_set(size()), v.get_mpi_communicator()); localized_vector.reinit(v, false, true); - Assert(localized_vector.size() == vec_size, - ExcDimensionMismatch(localized_vector.size(), vec_size)); + Assert(localized_vector.size() == size(), + ExcDimensionMismatch(localized_vector.size(), size())); // get a representation of the vector // and copy it @@ -941,7 +863,7 @@ Vector::operator=(const TrilinosWrappers::MPI::Vector &v) int ierr = localized_vector.trilinos_vector().ExtractView(&start_ptr); AssertThrow(ierr == 0, ExcTrilinosError(ierr)); - std::copy(start_ptr[0], start_ptr[0] + vec_size, begin()); + std::copy(start_ptr[0], start_ptr[0] + size(), begin()); } return *this; @@ -954,8 +876,7 @@ template bool Vector::operator==(const Vector &v) const { - Assert(vec_size != 0, ExcEmptyObject()); - Assert(vec_size == v.size(), ExcDimensionMismatch(vec_size, v.size())); + Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size())); // compare the two vector. we have to // convert the elements of v to the type of @@ -963,7 +884,7 @@ Vector::operator==(const Vector &v) const // because // operator==(complex,complex) // is not defined by default - for (size_type i = 0; i < vec_size; ++i) + for (size_type i = 0; i < size(); ++i) if (values[i] != Number(v.values[i])) return false; @@ -976,7 +897,7 @@ template void Vector::print(const char *format) const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); for (size_type j = 0; j < size(); ++j) internal::VectorOperations::print(values[j], format); @@ -992,7 +913,7 @@ Vector::print(std::ostream & out, const bool scientific, const bool across) const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); AssertThrow(out, ExcIO()); std::ios::fmtflags old_flags = out.flags(); @@ -1026,7 +947,7 @@ Vector::print(LogStream & out, const unsigned int width, const bool across) const { - Assert(vec_size != 0, ExcEmptyObject()); + Assert(size() != 0, ExcEmptyObject()); if (across) for (size_type i = 0; i < size(); ++i) @@ -1121,28 +1042,10 @@ template std::size_t Vector::memory_consumption() const { - return sizeof(*this) + (max_vec_size * sizeof(Number)); + return sizeof(*this) + values.memory_consumption() - sizeof(values); } - -template -void -Vector::allocate(const size_type copy_n_el) -{ - // allocate memory with the proper alignment requirements of 64 bytes - Number *new_values; - Utilities::System::posix_memalign(reinterpret_cast(&new_values), - 64, - sizeof(Number) * max_vec_size); - // copy: - for (size_type i = 0; i < copy_n_el; ++i) - new_values[i] = values[i]; - values.reset(new_values); -} - - - DEAL_II_NAMESPACE_CLOSE #endif diff --git a/tests/serialization/vector.output b/tests/serialization/vector.output index 2fc4d411cf..d51d6c2c02 100644 --- a/tests/serialization/vector.output +++ b/tests/serialization/vector.output @@ -1,6 +1,6 @@ -DEAL::0 0 0 0 5 5 0 1 2 3 4 +DEAL::0 0 0 0 0 0 5 0 1 2 3 4 -DEAL::0 0 0 0 5 5 0 1 2 3 4 +DEAL::0 0 0 0 0 0 5 0 1 2 3 4 DEAL::OK -- 2.39.5