From: Daniel Arndt Date: Thu, 16 Nov 2017 14:41:37 +0000 (+0100) Subject: Introduce ArrayView::data() and use where appropriate X-Git-Tag: v9.0.0-rc1~760^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3d4d3135556e877657a8909e3323be50405af7e4;p=dealii.git Introduce ArrayView::data() and use where appropriate --- diff --git a/include/deal.II/base/array_view.h b/include/deal.II/base/array_view.h index dbeb9569b6..a6e6110b1d 100644 --- a/include/deal.II/base/array_view.h +++ b/include/deal.II/base/array_view.h @@ -193,6 +193,12 @@ public: */ std::size_t size() const; + /** + * Return a pointer to the underlying array serving as element storage. + * In case the container is empty a nullptr is returned. + */ + value_type *data() const noexcept; + /** * Return an iterator pointing to the beginning of the array view. */ @@ -256,9 +262,8 @@ ArrayView::ArrayView() template inline -ArrayView:: -ArrayView(value_type *starting_element, - const std::size_t n_elements) +ArrayView::ArrayView(value_type *starting_element, + const std::size_t n_elements) : starting_element (starting_element), n_elements(n_elements) @@ -319,7 +324,7 @@ inline bool ArrayView::operator == (const ArrayView &other_view) const { - return (other_view.begin() == starting_element) + return (other_view.data() == starting_element) && (other_view.size() == n_elements); } @@ -331,7 +336,7 @@ bool ArrayView::operator == (const ArrayView::type> &other_view) const { - return (other_view.begin() == starting_element) + return (other_view.data() == starting_element) && (other_view.size() == n_elements); } @@ -347,6 +352,19 @@ ArrayView::operator != (const ArrayView &other_vi +template +inline +typename ArrayView::value_type * +ArrayView::data() const noexcept +{ + if (n_elements==0) + return nullptr; + else + return starting_element; +} + + + template inline bool @@ -375,6 +393,7 @@ ArrayView::begin() const } + template inline typename ArrayView::iterator @@ -383,6 +402,8 @@ ArrayView::end() const return starting_element + n_elements; } + + template inline typename ArrayView::const_iterator @@ -392,6 +413,7 @@ ArrayView::cbegin() const } + template inline typename ArrayView::const_iterator @@ -401,6 +423,7 @@ ArrayView::cend() const } + template inline typename ArrayView::value_type & diff --git a/include/deal.II/base/mpi.templates.h b/include/deal.II/base/mpi.templates.h index 4912534261..4ad8353035 100644 --- a/include/deal.II/base/mpi.templates.h +++ b/include/deal.II/base/mpi.templates.h @@ -127,12 +127,12 @@ namespace Utilities // implementations of MPI-2. It is not needed as // of MPI-3 and we should remove it at some // point in the future. - const_cast(static_cast(values.begin())) + const_cast(static_cast(values.data())) : MPI_IN_PLACE, - static_cast(output.begin()), + static_cast(output.data()), static_cast(values.size()), - internal::mpi_type_id(values.begin()), + internal::mpi_type_id(values.data()), mpi_op, mpi_communicator); AssertThrowMPI(ierr); @@ -167,10 +167,10 @@ namespace Utilities // implementations of MPI-2. It is not needed as // of MPI-3 and we should remove it at some // point in the future. - const_cast(static_cast(values.begin())) + const_cast(static_cast(values.data())) : MPI_IN_PLACE, - static_cast(output.begin()), + static_cast(output.data()), static_cast(values.size()*2), internal::mpi_type_id(static_cast(nullptr)), mpi_op, diff --git a/include/deal.II/base/partitioner.templates.h b/include/deal.II/base/partitioner.templates.h index 398a8bb874..5ba72ce87d 100644 --- a/include/deal.II/base/partitioner.templates.h +++ b/include/deal.II/base/partitioner.templates.h @@ -69,9 +69,9 @@ namespace Utilities const bool use_larger_set = (n_ghost_indices_in_larger_set > n_ghost_indices() && ghost_array.size() == n_ghost_indices_in_larger_set); Number *ghost_array_ptr = use_larger_set ? - ghost_array.begin()+ + ghost_array.data()+ n_ghost_indices_in_larger_set-n_ghost_indices() - : ghost_array.begin(); + : ghost_array.data(); for (unsigned int i=0; i(import_targets_data[i].second)* @@ -238,7 +238,7 @@ namespace Utilities // in case we want to import only from a subset of the ghosts we want to // move the data to send to the front of the array AssertIndexRange(n_ghost_indices(), n_ghost_indices_in_larger_set+1); - Number *ghost_array_ptr = ghost_array.begin(); + Number *ghost_array_ptr = ghost_array.data(); for (unsigned int i=0; ifirst) + if (ghost_array_ptr + offset != ghost_array.data() + my_ghosts->first) for (unsigned int j=my_ghosts->first; jsecond; ++j, ++offset) { ghost_array_ptr[offset] = ghost_array[j]; @@ -334,7 +334,7 @@ namespace Utilities "vector_operation argument was passed to " "import_from_ghosted_array_start as is passed " "to import_from_ghosted_array_finish.")); - std::memset(ghost_array.begin(), 0, sizeof(Number)*ghost_array.size()); + std::memset(ghost_array.data(), 0, sizeof(Number)*ghost_array.size()); return; } #endif @@ -357,7 +357,7 @@ namespace Utilities MPI_STATUSES_IGNORE); AssertThrowMPI(ierr); - const Number *read_position = temporary_storage.begin(); + const Number *read_position = temporary_storage.data(); std::vector >::const_iterator my_imports = import_indices_data.begin(); @@ -386,7 +386,7 @@ namespace Utilities typename LinearAlgebra::distributed::Vector:: ExcNonMatchingElements(*read_position, locally_owned_array[j], my_pid)); - AssertDimension(read_position-temporary_storage.begin(), n_import_indices()); + AssertDimension(read_position-temporary_storage.data(), n_import_indices()); } // wait for the send operations to complete @@ -405,7 +405,7 @@ namespace Utilities if (ghost_array.size()>0) { Assert(ghost_array.begin()!=nullptr, ExcInternalError()); - std::memset(ghost_array.begin(), 0, sizeof(Number)*n_ghost_indices()); + std::memset(ghost_array.data(), 0, sizeof(Number)*n_ghost_indices()); } // clear the compress requests diff --git a/include/deal.II/lac/tensor_product_matrix.h b/include/deal.II/lac/tensor_product_matrix.h index 5fbad76fa4..b60716e7b1 100644 --- a/include/deal.II/lac/tensor_product_matrix.h +++ b/include/deal.II/lac/tensor_product_matrix.h @@ -520,7 +520,7 @@ TensorProductMatrixSymmetricSumBase eval(AlignedVector(), AlignedVector(), AlignedVector(), mass_matrix[0].n_rows()-1, mass_matrix[0].n_rows()); Number *t = tmp_array.begin(); - const Number *src = src_view.begin(); + const Number *src = src_view.data(); Number *dst = &(dst_view[0]); // NOTE: dof_to_quad has to be interpreted as 'dof to eigenvalue index' diff --git a/source/grid/manifold.cc b/source/grid/manifold.cc index 67c8c4ce81..687eba0971 100644 --- a/source/grid/manifold.cc +++ b/source/grid/manifold.cc @@ -571,7 +571,7 @@ get_new_points (const ArrayView> &surrounding_points, // check whether periodicity shifts some of the points. Only do this if // necessary to avoid memory allocation - const Point *surrounding_points_start = surrounding_points.begin(); + const Point *surrounding_points_start = surrounding_points.data(); boost::container::small_vector, 200> modified_points; bool adjust_periodicity = false;