]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Replace C-style arrays in MatrixFree 10895/head
authorPeter Munch <peterrmuench@gmail.com>
Sat, 5 Sep 2020 19:02:32 +0000 (21:02 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Sat, 5 Sep 2020 19:02:32 +0000 (21:02 +0200)
include/deal.II/matrix_free/dof_info.h
include/deal.II/matrix_free/dof_info.templates.h
include/deal.II/matrix_free/evaluation_kernels.h
include/deal.II/matrix_free/face_info.h
include/deal.II/matrix_free/mapping_info.h
include/deal.II/matrix_free/shape_info.h

index cef97730b9fb47b98c745468513a26f0393b0ccd..36002f9a90b7d97f7b6ae32c3be2b15531759d29 100644 (file)
@@ -433,7 +433,7 @@ namespace internal
        * as interior (0), the faces decorated with as exterior (1), and the
        * cells (2) according to CellOrFaceAccess.
        */
-      std::vector<IndexStorageVariants> index_storage_variants[3];
+      std::array<std::vector<IndexStorageVariants>, 3> index_storage_variants;
 
       /**
        * Stores the rowstart indices of the compressed row storage in the @p
@@ -485,7 +485,7 @@ namespace internal
        * as interior (0), the faces decorated with as exterior (1), and the
        * cells (2) according to CellOrFaceAccess.
        */
-      std::vector<unsigned int> dof_indices_contiguous[3];
+      std::array<std::vector<unsigned int>, 3> dof_indices_contiguous;
 
       /**
        * Compressed index storage for faster access than through @p
@@ -495,7 +495,7 @@ namespace internal
        * as minus (0), the faces decorated with as plus (1), and the cells
        * (2).
        */
-      std::vector<unsigned int> dof_indices_interleave_strides[3];
+      std::array<std::vector<unsigned int>, 3> dof_indices_interleave_strides;
 
       /**
        * Caches the number of indices filled when vectorizing. This
@@ -506,7 +506,7 @@ namespace internal
        * as interior (0), the faces decorated with as exterior (1), and the
        * cells (2) according to CellOrFaceAccess.
        */
-      std::vector<unsigned char> n_vectorization_lanes_filled[3];
+      std::array<std::vector<unsigned char>, 3> n_vectorization_lanes_filled;
 
       /**
        * This stores the parallel partitioning that can be used to set up
index 75761d73db836cb766172c63f14cf6482b303298..f8d115f6d89406ae08c5325ccb666605cad395e3 100644 (file)
@@ -1066,7 +1066,8 @@ namespace internal
       for (unsigned int face = 0; face < faces.size(); ++face)
         {
           auto face_computation = [&](const DoFAccessIndex face_index,
-                                      const unsigned int * cell_indices_face) {
+                                      const std::array<unsigned int, length>
+                                        &cell_indices_face) {
             bool is_contiguous      = false;
             bool is_interleaved     = false;
             bool needs_full_storage = false;
index 1c9121b84779e94e3667d83026041446ca089438..0ffed895902abb7dc6385a86934749dbbe17f1ec 100644 (file)
@@ -1784,14 +1784,14 @@ namespace internal
   private:
     template <bool do_evaluate, bool add_into_output, int face_direction = 0>
     static void
-    interpolate_generic(const Number *               input,
-                        Number *                     output,
-                        const bool                   do_gradients,
-                        const unsigned int           face_no,
-                        const unsigned int           n_points_1d,
-                        const AlignedVector<Number> *shape_data,
-                        const unsigned int           dofs_per_component_on_cell,
-                        const unsigned int           dofs_per_component_on_face)
+    interpolate_generic(const Number *     input,
+                        Number *           output,
+                        const bool         do_gradients,
+                        const unsigned int face_no,
+                        const unsigned int n_points_1d,
+                        const std::array<AlignedVector<Number>, 2> &shape_data,
+                        const unsigned int dofs_per_component_on_cell,
+                        const unsigned int dofs_per_component_on_face)
     {
       if (face_direction == face_no / 2)
         {
index 43d8d820e822f8f2bab3cf0c9430255e0c09d650..b411c27d7e1bfea5a66c7b5dcbb8892bb3255c18 100644 (file)
@@ -58,7 +58,7 @@ namespace internal
        * numbers of the cells on the logical "interior" side of the face which
        * is aligned to the direction of FEEvaluation::get_normal_vector().
        */
-      unsigned int cells_interior[vectorization_width];
+      std::array<unsigned int, vectorization_width> cells_interior;
 
       /**
        * Indices of the faces in the current face batch as compared to the
@@ -73,7 +73,7 @@ namespace internal
        * For boundary faces, the numbers are set to
        * `numbers::invalid_unsigned_int`.
        */
-      unsigned int cells_exterior[vectorization_width];
+      std::array<unsigned int, vectorization_width> cells_exterior;
 
       /**
        * Index of the face between 0 and GeometryInfo::faces_per_cell within
index 4d7c9f64c5814fb409d76e61cb9a83fe6829ff70..b2ba970de7c856701cea70eb23f9592401192267 100644 (file)
@@ -220,7 +220,8 @@ namespace internal
        * but the default case (cell integrals or boundary integrals) only
        * fills the zeroth component and ignores the first one.
        */
-      AlignedVector<Tensor<2, spacedim, VectorizedArrayType>> jacobians[2];
+      std::array<AlignedVector<Tensor<2, spacedim, VectorizedArrayType>>, 2>
+        jacobians;
 
       /**
        * The storage of the gradients of the inverse Jacobian
@@ -237,10 +238,12 @@ namespace internal
        * but the default case (cell integrals or boundary integrals) only
        * fills the zeroth component and ignores the first one.
        */
-      AlignedVector<Tensor<1,
-                           spacedim *(spacedim + 1) / 2,
-                           Tensor<1, spacedim, VectorizedArrayType>>>
-        jacobian_gradients[2];
+      std::array<
+        AlignedVector<Tensor<1,
+                             spacedim *(spacedim + 1) / 2,
+                             Tensor<1, spacedim, VectorizedArrayType>>>,
+        2>
+        jacobian_gradients;
 
       /**
        * Stores the Jacobian transformations times the normal vector (this
@@ -249,8 +252,8 @@ namespace internal
        *
        * Indexed by @p data_index_offsets.
        */
-      AlignedVector<Tensor<1, spacedim, VectorizedArrayType>>
-        normals_times_jacobians[2];
+      std::array<AlignedVector<Tensor<1, spacedim, VectorizedArrayType>>, 2>
+        normals_times_jacobians;
 
       /**
        * Stores the index offset of a particular cell into the quadrature
index b30412de0862f797be8ff99933a8a20a1c78d7e5..0637a5e775be7450e0cbf8f29123f763f860b757 100644 (file)
@@ -227,7 +227,7 @@ namespace internal
        * (the vertices) in one data structure. Sorting is first the values,
        * then gradients, then second derivatives.
        */
-      AlignedVector<Number> shape_data_on_face[2];
+      std::array<AlignedVector<Number>, 2> shape_data_on_face;
 
       /**
        * Collects all data of 1D nodal shape values (defined by the Lagrange
@@ -239,25 +239,25 @@ namespace internal
        *
        * @note In contrast to shape_data_on_face, only the vales are evaluated.
        */
-      AlignedVector<Number> quadrature_data_on_face[2];
+      std::array<AlignedVector<Number>, 2> quadrature_data_on_face;
 
       /**
        * Stores one-dimensional values of shape functions on subface. Since
        * there are two subfaces, store two variants.
        */
-      AlignedVector<Number> values_within_subface[2];
+      std::array<AlignedVector<Number>, 2> values_within_subface;
 
       /**
        * Stores one-dimensional gradients of shape functions on subface. Since
        * there are two subfaces, store two variants.
        */
-      AlignedVector<Number> gradients_within_subface[2];
+      std::array<AlignedVector<Number>, 2> gradients_within_subface;
 
       /**
        * Stores one-dimensional gradients of shape functions on subface. Since
        * there are two subfaces, store two variants.
        */
-      AlignedVector<Number> hessians_within_subface[2];
+      std::array<AlignedVector<Number>, 2> hessians_within_subface;
 
       /**
        * We store a copy of the one-dimensional quadrature formula

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.