]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Introduce ArrayView::data() and use where appropriate
authorDaniel Arndt <d.arndt@math.uni-goettingen.de>
Thu, 16 Nov 2017 14:41:37 +0000 (15:41 +0100)
committerDaniel Arndt <d.arndt@math.uni-goettingen.de>
Fri, 17 Nov 2017 01:22:09 +0000 (02:22 +0100)
include/deal.II/base/array_view.h
include/deal.II/base/mpi.templates.h
include/deal.II/base/partitioner.templates.h
include/deal.II/lac/tensor_product_matrix.h
source/grid/manifold.cc

index dbeb9569b659485476afd3e0dd5d965c8d3375d5..a6e6110b1dc3f7a46989c9568bc322af0fd92a45 100644 (file)
@@ -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<ElementType>::ArrayView()
 
 template <typename ElementType>
 inline
-ArrayView<ElementType>::
-ArrayView(value_type        *starting_element,
-          const std::size_t  n_elements)
+ArrayView<ElementType>::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<ElementType>::operator == (const ArrayView<const value_type> &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<ElementType>::operator ==
 (const ArrayView<typename std::remove_cv<value_type>::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<ElementType>::operator != (const ArrayView<const value_type> &other_vi
 
 
 
+template <typename ElementType>
+inline
+typename ArrayView<ElementType>::value_type *
+ArrayView<ElementType>::data() const noexcept
+{
+  if (n_elements==0)
+    return nullptr;
+  else
+    return starting_element;
+}
+
+
+
 template <typename ElementType>
 inline
 bool
@@ -375,6 +393,7 @@ ArrayView<ElementType>::begin() const
 }
 
 
+
 template <typename ElementType>
 inline
 typename ArrayView<ElementType>::iterator
@@ -383,6 +402,8 @@ ArrayView<ElementType>::end() const
   return starting_element + n_elements;
 }
 
+
+
 template <typename ElementType>
 inline
 typename ArrayView<ElementType>::const_iterator
@@ -392,6 +413,7 @@ ArrayView<ElementType>::cbegin() const
 }
 
 
+
 template <typename ElementType>
 inline
 typename ArrayView<ElementType>::const_iterator
@@ -401,6 +423,7 @@ ArrayView<ElementType>::cend() const
 }
 
 
+
 template <typename ElementType>
 inline
 typename ArrayView<ElementType>::value_type &
index 4912534261e5143f9586630bdd06a165246cd3fc..4ad83530356a0410c44662ec98572c1f48b6029c 100644 (file)
@@ -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<void *>(static_cast<const void *>(values.begin()))
+                              const_cast<void *>(static_cast<const void *>(values.data()))
                               :
                               MPI_IN_PLACE,
-                              static_cast<void *>(output.begin()),
+                              static_cast<void *>(output.data()),
                               static_cast<int>(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<void *>(static_cast<const void *>(values.begin()))
+                              const_cast<void *>(static_cast<const void *>(values.data()))
                               :
                               MPI_IN_PLACE,
-                              static_cast<void *>(output.begin()),
+                              static_cast<void *>(output.data()),
                               static_cast<int>(values.size()*2),
                               internal::mpi_type_id(static_cast<T *>(nullptr)),
                               mpi_op,
index 398a8bb87457c3513d45145bb21a1c46442fe57e..5ba72ce87dfba922a5462d5f88c823318742d856 100644 (file)
@@ -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<n_ghost_targets; i++)
         {
@@ -88,7 +88,7 @@ namespace Utilities
           ghost_array_ptr += ghost_targets()[i].second;
         }
 
-      Number *temp_array_ptr = temporary_storage.begin();
+      Number *temp_array_ptr = temporary_storage.data();
       for (unsigned int i=0; i<n_import_targets; i++)
         {
           // copy the data to be sent to the import_data field
@@ -213,7 +213,7 @@ namespace Utilities
       requests.resize (n_import_targets + n_ghost_targets);
 
       // initiate the receive operations
-      Number *temp_array_ptr = temporary_storage.begin();
+      Number *temp_array_ptr = temporary_storage.data();
       for (unsigned int i=0; i<n_import_targets; i++)
         {
           AssertThrow (static_cast<std::size_t>(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; i<n_ghost_targets; i++)
         {
           // in case we only sent a subset of indices, we now need to move the data
@@ -251,7 +251,7 @@ namespace Utilities
               end_my_ghosts = ghost_indices_subset_data.begin()+ghost_indices_subset_chunks_by_rank_data[i+1];
               unsigned int offset = 0;
               for ( ; my_ghosts != end_my_ghosts; ++my_ghosts)
-                if (ghost_array_ptr + offset != ghost_array.begin() + my_ghosts->first)
+                if (ghost_array_ptr + offset != ghost_array.data() + my_ghosts->first)
                   for (unsigned int j=my_ghosts->first; j<my_ghosts->second; ++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<std::pair<unsigned int, unsigned int> >::const_iterator
           my_imports = import_indices_data.begin();
 
@@ -386,7 +386,7 @@ namespace Utilities
                        typename LinearAlgebra::distributed::Vector<Number>::
                        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
index 5fbad76fa4c96023d5b56075cb074257a4041cb7..b60716e7b157cbe9cc319910b4211f46c8ef0be6 100644 (file)
@@ -520,7 +520,7 @@ TensorProductMatrixSymmetricSumBase<dim,Number,size>
   eval(AlignedVector<Number>(), AlignedVector<Number>(),
        AlignedVector<Number>(), 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'
index 67c8c4ce81e0a7a6b9f1744777825f3b4b9a748d..687eba0971037516baeb6249d19d8c4875cde1b4 100644 (file)
@@ -571,7 +571,7 @@ get_new_points (const ArrayView<const Point<spacedim>> &surrounding_points,
 
   // check whether periodicity shifts some of the points. Only do this if
   // necessary to avoid memory allocation
-  const Point<spacedim> *surrounding_points_start = surrounding_points.begin();
+  const Point<spacedim> *surrounding_points_start = surrounding_points.data();
 
   boost::container::small_vector<Point<spacedim>, 200> modified_points;
   bool adjust_periodicity = false;

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.