From 7ad068984b29f5712dd2bed5742459a9e2030ac3 Mon Sep 17 00:00:00 2001 From: David Wells Date: Sat, 7 Jan 2023 19:17:36 -0500 Subject: [PATCH] Avoid using Tensor::begin_raw() et al. These are deprecated and in several cases we have better choices anyway. --- include/deal.II/base/bounding_box.h | 11 +++++----- source/base/quadrature_lib.cc | 17 +++++++-------- source/base/tensor.cc | 17 ++++++++++++--- source/fe/fe_simplex_p.cc | 7 ++----- source/fe/fe_simplex_p_bubbles.cc | 8 +++---- tests/base/quadrature_point_data_04.cc | 6 ++++-- .../symengine_tensor_operations_04.cc | 21 ++++++++++++------- 7 files changed, 49 insertions(+), 38 deletions(-) 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::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::infinity()); - std::fill(max.begin_raw(), - max.end_raw(), - -std::numeric_limits::infinity()); + for (unsigned int d = 0; d < spacedim; ++d) + { + min[d] = std::numeric_limits::infinity(); + max[d] = -std::numeric_limits::infinity(); + } for (const Point &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::QWitherdenVincentSimplex( const double volume = (dim == 2 ? 1.0 / 2.0 : 1.0 / 6.0); this->weights.emplace_back(volume * b_weights[permutation_n]); Point 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 #include #include @@ -45,22 +46,32 @@ namespace const types::blas_int lwork = 5 * dim; std::array work; types::blas_int info; + constexpr std::size_t size = + Tensor<2, dim, Number>::n_independent_components; + std::array A_array; + A_in_VT_out.unroll(A_array.begin(), A_array.end()); + std::array 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 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().template barycenter()); 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 fe_p(degree); std::vector> points = fe_p.get_unit_support_points(); - Point centroid; - std::fill(centroid.begin_raw(), centroid.end_raw(), 1.0 / double(dim + 1)); + const Point centroid = + fe_p.reference_cell().template barycenter(); switch (dim) { @@ -117,8 +117,8 @@ namespace FE_P_BubblesImplementation BarycentricPolynomials get_basis(const unsigned int degree) { - Point centroid; - std::fill(centroid.begin_raw(), centroid.end_raw(), 1.0 / double(dim + 1)); + const Point centroid = + ReferenceCells::get_simplex().template barycenter(); auto M = [](const unsigned int d) { return BarycentricPolynomial::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 &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 &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; 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("a"); @@ -128,10 +130,13 @@ test_symmetric_tensor() using Tensor_t = SymmetricTensor; 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("a"); -- 2.39.5