From: David Wells <drwells@email.unc.edu>
Date: Sun, 8 Jan 2023 00:17:36 +0000 (-0500)
Subject: Avoid using Tensor::begin_raw() et al.
X-Git-Tag: v9.5.0-rc1~674^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F14647%2Fhead;p=dealii.git

Avoid using Tensor::begin_raw() et al.

These are deprecated and in several cases we have better choices anyway.
---

diff --git a/include/deal.II/base/bounding_box.h b/include/deal.II/base/bounding_box.h
index e2bc943112..6ea5bc3354 100644
--- a/include/deal.II/base/bounding_box.h
+++ b/include/deal.II/base/bounding_box.h
@@ -432,12 +432,11 @@ inline BoundingBox<spacedim, Number>::BoundingBox(const Container &points)
     {
       auto &min = boundary_points.first;
       auto &max = boundary_points.second;
-      std::fill(min.begin_raw(),
-                min.end_raw(),
-                std::numeric_limits<Number>::infinity());
-      std::fill(max.begin_raw(),
-                max.end_raw(),
-                -std::numeric_limits<Number>::infinity());
+      for (unsigned int d = 0; d < spacedim; ++d)
+        {
+          min[d] = std::numeric_limits<Number>::infinity();
+          max[d] = -std::numeric_limits<Number>::infinity();
+        }
 
       for (const Point<spacedim, Number> &point : points)
         for (unsigned int d = 0; d < spacedim; ++d)
diff --git a/source/base/quadrature_lib.cc b/source/base/quadrature_lib.cc
index 7803cb6196..df8cb2efac 100644
--- a/source/base/quadrature_lib.cc
+++ b/source/base/quadrature_lib.cc
@@ -615,13 +615,11 @@ template <>
 unsigned int
 QGaussOneOverR<2>::quad_size(const Point<2> &singularity, const unsigned int n)
 {
-  const double eps = 1e-8;
-  const bool   on_edge =
-    std::any_of(singularity.begin_raw(),
-                singularity.end_raw(),
-                [eps](double coord) {
-                  return std::abs(coord) < eps || std::abs(coord - 1.) < eps;
-                });
+  const double eps     = 1e-8;
+  bool         on_edge = false;
+  for (unsigned int d = 0; d < 2; ++d)
+    on_edge = on_edge || (std::abs(singularity[d]) < eps ||
+                          std::abs(singularity[d] - 1.0) < eps);
   const bool on_vertex =
     on_edge &&
     std::abs((singularity - Point<2>(.5, .5)).norm_square() - .5) < eps;
@@ -2022,9 +2020,8 @@ QWitherdenVincentSimplex<dim>::QWitherdenVincentSimplex(
           const double volume = (dim == 2 ? 1.0 / 2.0 : 1.0 / 6.0);
           this->weights.emplace_back(volume * b_weights[permutation_n]);
           Point<dim> c_point;
-          std::copy(b_point.begin(),
-                    b_point.begin() + dim,
-                    c_point.begin_raw());
+          for (int d = 0; d < dim; ++d)
+            c_point[d] = b_point[d];
           this->quadrature_points.emplace_back(c_point);
         }
     }
diff --git a/source/base/tensor.cc b/source/base/tensor.cc
index ad70671b68..25cec533ef 100644
--- a/source/base/tensor.cc
+++ b/source/base/tensor.cc
@@ -13,6 +13,7 @@
 //
 // ---------------------------------------------------------------------
 
+#include <deal.II/base/array_view.h>
 #include <deal.II/base/tensor.h>
 
 #include <deal.II/lac/exceptions.h>
@@ -45,22 +46,32 @@ namespace
     const types::blas_int     lwork = 5 * dim;
     std::array<Number, lwork> work;
     types::blas_int           info;
+    constexpr std::size_t     size =
+      Tensor<2, dim, Number>::n_independent_components;
+    std::array<Number, size> A_array;
+    A_in_VT_out.unroll(A_array.begin(), A_array.end());
+    std::array<Number, size> U_array;
+    U.unroll(U_array.begin(), U_array.end());
     gesvd(&LAPACKSupport::O, // replace VT in place
           &LAPACKSupport::A,
           &N,
           &N,
-          A_in_VT_out.begin_raw(),
+          A_array.data(),
           &N,
           S.data(),
-          A_in_VT_out.begin_raw(),
+          A_array.data(),
           &N,
-          U.begin_raw(),
+          U_array.data(),
           &N,
           work.data(),
           &lwork,
           &info);
     Assert(info == 0, LAPACKSupport::ExcErrorCode("gesvd", info));
     Assert(S.back() / S.front() > 1.e-10, LACExceptions::ExcSingular());
+
+    A_in_VT_out =
+      Tensor<2, dim, Number>(make_array_view(A_array.begin(), A_array.end()));
+    U = Tensor<2, dim, Number>(make_array_view(U_array.begin(), U_array.end()));
   }
 } // namespace
 
diff --git a/source/fe/fe_simplex_p.cc b/source/fe/fe_simplex_p.cc
index 155cf32169..18ab94784d 100644
--- a/source/fe/fe_simplex_p.cc
+++ b/source/fe/fe_simplex_p.cc
@@ -95,11 +95,8 @@ namespace
     // centroid and only the centroid
     if (degree == 0)
       {
-        Point<dim> centroid;
-        std::fill(centroid.begin_raw(),
-                  centroid.end_raw(),
-                  1.0 / double(dim + 1));
-        unit_points.emplace_back(centroid);
+        unit_points.emplace_back(
+          ReferenceCells::get_simplex<dim>().template barycenter<dim>());
         return unit_points;
       }
 
diff --git a/source/fe/fe_simplex_p_bubbles.cc b/source/fe/fe_simplex_p_bubbles.cc
index 990c59cd59..260e5a3847 100644
--- a/source/fe/fe_simplex_p_bubbles.cc
+++ b/source/fe/fe_simplex_p_bubbles.cc
@@ -67,8 +67,8 @@ namespace FE_P_BubblesImplementation
     FE_SimplexP<dim>        fe_p(degree);
     std::vector<Point<dim>> points = fe_p.get_unit_support_points();
 
-    Point<dim> centroid;
-    std::fill(centroid.begin_raw(), centroid.end_raw(), 1.0 / double(dim + 1));
+    const Point<dim> centroid =
+      fe_p.reference_cell().template barycenter<dim>();
 
     switch (dim)
       {
@@ -117,8 +117,8 @@ namespace FE_P_BubblesImplementation
   BarycentricPolynomials<dim>
   get_basis(const unsigned int degree)
   {
-    Point<dim> centroid;
-    std::fill(centroid.begin_raw(), centroid.end_raw(), 1.0 / double(dim + 1));
+    const Point<dim> centroid =
+      ReferenceCells::get_simplex<dim>().template barycenter<dim>();
 
     auto M = [](const unsigned int d) {
       return BarycentricPolynomial<dim, double>::monomial(d);
diff --git a/tests/base/quadrature_point_data_04.cc b/tests/base/quadrature_point_data_04.cc
index 3f00baad22..c7a7d8b7e4 100644
--- a/tests/base/quadrature_point_data_04.cc
+++ b/tests/base/quadrature_point_data_04.cc
@@ -66,13 +66,15 @@ struct Mat1 : MaterialBase
   pack_values(std::vector<double> &values) const final
   {
     AssertDimension(values.size(), number_of_values());
-    std::copy(pt.begin_raw(), pt.end_raw(), values.begin());
+    for (unsigned int d = 0; d < 2; ++d)
+      values[d] = pt[d];
   }
   virtual void
   unpack_values(const std::vector<double> &values) final
   {
     AssertDimension(values.size(), number_of_values());
-    std::copy(values.cbegin(), values.cend(), pt.begin_raw());
+    for (unsigned int d = 0; d < 2; ++d)
+      pt[d] = values[d];
   }
 };
 
diff --git a/tests/symengine/symengine_tensor_operations_04.cc b/tests/symengine/symengine_tensor_operations_04.cc
index 06af8c9a2f..71510f7676 100644
--- a/tests/symengine/symengine_tensor_operations_04.cc
+++ b/tests/symengine/symengine_tensor_operations_04.cc
@@ -104,10 +104,12 @@ test_tensor()
   using Tensor_t           = Tensor<rank, dim, double>;
 
   Tensor_t t_a, t_b;
-  for (auto it = t_a.begin_raw(); it != t_a.end_raw(); ++it)
-    *it = 1.0;
-  for (auto it = t_b.begin_raw(); it != t_b.end_raw(); ++it)
-    *it = 2.0;
+  for (unsigned int i = 0; i < Tensor_t::n_independent_components; ++i)
+    {
+      const auto index = t_a.unrolled_to_component_index(i);
+      t_a[index]       = 1.0;
+      t_b[index]       = 2.0;
+    }
 
   const Tensor_SD_number_t symb_t_a =
     SD::make_tensor_of_symbols<rank, dim>("a");
@@ -128,10 +130,13 @@ test_symmetric_tensor()
   using Tensor_t           = SymmetricTensor<rank, dim, double>;
 
   Tensor_t t_a, t_b;
-  for (auto it = t_a.begin_raw(); it != t_a.end_raw(); ++it)
-    *it = 1.0;
-  for (auto it = t_b.begin_raw(); it != t_b.end_raw(); ++it)
-    *it = 2.0;
+  for (unsigned int i = 0; i < Tensor_t::n_independent_components; ++i)
+    {
+      const auto index = t_a.unrolled_to_component_index(i);
+      t_a[index]       = 1.0;
+      t_b[index]       = 2.0;
+    }
+
 
   const Tensor_SD_number_t symb_t_a =
     SD::make_symmetric_tensor_of_symbols<rank, dim>("a");