]> https://gitweb.dealii.org/ - dealii.git/commitdiff
MatrixFree: Switch ShapeInfo from VectorizedArray<Number> to Number type 16985/head
authorMartin Kronbichler <martin.kronbichler@rub.de>
Sun, 5 May 2024 06:03:17 +0000 (08:03 +0200)
committerMartin Kronbichler <martin.kronbichler@rub.de>
Thu, 9 May 2024 10:12:34 +0000 (12:12 +0200)
15 files changed:
examples/step-76/step-76.cc
include/deal.II/fe/mapping_q.h
include/deal.II/matrix_free/constraint_info.h
include/deal.II/matrix_free/evaluation_kernels.h
include/deal.II/matrix_free/evaluation_kernels_face.h
include/deal.II/matrix_free/evaluation_kernels_hanging_nodes.h
include/deal.II/matrix_free/evaluation_template_factory.h
include/deal.II/matrix_free/evaluation_template_factory_hanging_nodes.templates.h
include/deal.II/matrix_free/fe_evaluation.h
include/deal.II/matrix_free/fe_evaluation_data.h
include/deal.II/matrix_free/mapping_info.templates.h
include/deal.II/matrix_free/matrix_free.h
tests/matrix_free/fe_evaluation_shift.cc
tests/matrix_free/hanging_node_kernels_01.cc
tests/performance/timing_navier_stokes.cc

index f01c40479d7426685ed151268637abf00b07607f..90f67c66e019011629701fa49a09945187d32db5 100644 (file)
@@ -663,10 +663,11 @@ namespace Euler_DG
           dim,
           n_points_1d,
           n_points_1d,
-          VectorizedArrayType>
-          eval(AlignedVector<VectorizedArrayType>(),
+          VectorizedArrayType,
+          Number>
+          eval({},
                data.get_shape_info().data[0].shape_gradients_collocation_eo,
-               AlignedVector<VectorizedArrayType>());
+               {});
 
         AlignedVector<VectorizedArrayType> buffer(phi.static_n_q_points *
                                                   phi.n_components);
index a440434ae6186998fbf30bdb6b04b3dc875c9b04..69b77a776e88ceef371ce40782ea525e00f6cbd8 100644 (file)
@@ -405,7 +405,7 @@ public:
      * we need to store the evaluations of the 1d polynomials at
      * the 1d quadrature points. That is what this variable is for.
      */
-    internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> shape_info;
+    internal::MatrixFreeFunctions::ShapeInfo<double> shape_info;
 
     /**
      * In case the quadrature rule given represents a tensor product
index 6348f96eaea2c95683432f3afeb228dac5904125..3dde53aabd999973f88ed4092ab6967448949ea1 100644 (file)
@@ -187,7 +187,7 @@ namespace internal
       std::vector<typename Number::value_type> constraint_pool_data;
       std::vector<unsigned int>                constraint_pool_row_index;
 
-      std::vector<ShapeInfo<Number>>          shape_infos;
+      std::vector<ShapeInfo<typename Number::value_type>> shape_infos;
       std::vector<compressed_constraint_kind> hanging_node_constraint_masks;
       std::vector<unsigned int>               active_fe_indices;
 
index 25419d5ded1100fe84626e181d80126affc781bf..3ec800f16ff26155d8afaac5d31ef44273c76bcb 100644 (file)
@@ -1285,7 +1285,7 @@ namespace internal
     // might have non-symmetric quadrature formula, so use the more
     // conservative 'evaluate_general' scheme rather than 'even_odd' as the
     // Hessians are not used very often
-    const MatrixFreeFunctions::UnivariateShapeData<Number> &data =
+    const MatrixFreeFunctions::UnivariateShapeData<Number2> &data =
       fe_eval.get_shape_info().data[0];
     AssertDimension(data.shape_gradients_collocation.size(),
                     data.n_q_points_1d * data.n_q_points_1d);
@@ -1356,7 +1356,7 @@ namespace internal
     using Number2 =
       typename FEEvaluationData<dim, Number, false>::shape_info_number_type;
 
-    const MatrixFreeFunctions::UnivariateShapeData<Number> &data =
+    const MatrixFreeFunctions::UnivariateShapeData<Number2> &data =
       fe_eval.get_shape_info().data[0];
     AssertDimension(data.shape_gradients_collocation.size(),
                     data.n_q_points_1d * data.n_q_points_1d);
index 974d273bbc31a0a0b3f16914601cd9305ce8687a..f6468be33cc07cfe83d742b211fe153c9c395038 100644 (file)
@@ -2429,10 +2429,10 @@ namespace internal
                 const auto face_no = fe_eval.get_face_no(v);
 
                 grad_weight[v] =
-                  shape_data.shape_data_on_face[0][fe_degree +
-                                                   (integrate ?
-                                                      (2 - (face_no % 2)) :
-                                                      (1 + (face_no % 2)))][0];
+                  shape_data
+                    .shape_data_on_face[0][fe_degree + (integrate ?
+                                                          (2 - (face_no % 2)) :
+                                                          (1 + (face_no % 2)))];
 
                 index_array_hermite[v] =
                   &fe_eval.get_shape_info().face_to_cell_index_hermite(face_no,
index 15959f99919970cadab3cfb48445c1c28cd63d7e..bbd49888b1d7805ceab0f4ec9fa76cf3cdfbc9a1 100644 (file)
@@ -88,15 +88,18 @@ namespace internal
     Number>
   {
   private:
-    template <int structdim, unsigned int direction, bool transpose>
+    template <int          structdim,
+              unsigned int direction,
+              bool         transpose,
+              typename Number2>
     static void
-    interpolate(const unsigned int             offset,
-                const unsigned int             outer_stride,
-                const unsigned int             given_degree,
-                const Number                   mask_weight,
-                const Number                   mask_write,
-                const Number *DEAL_II_RESTRICT weights,
-                Number *DEAL_II_RESTRICT       values)
+    interpolate(const unsigned int              offset,
+                const unsigned int              outer_stride,
+                const unsigned int              given_degree,
+                const Number                    mask_weight,
+                const Number                    mask_write,
+                const Number2 *DEAL_II_RESTRICT weights,
+                Number *DEAL_II_RESTRICT        values)
     {
       static constexpr unsigned int max_n_points_1D = 40;
 
@@ -178,19 +181,19 @@ namespace internal
     }
 
   public:
-    template <bool transpose>
+    template <bool transpose, typename Number2>
     static void
     run_internal(
-      const unsigned int                            n_components,
-      const MatrixFreeFunctions::ShapeInfo<Number> &shape_info,
+      const unsigned int                             n_components,
+      const MatrixFreeFunctions::ShapeInfo<Number2> &shape_info,
       const std::array<MatrixFreeFunctions::compressed_constraint_kind,
-                       Number::size()>             &constraint_mask,
-      Number                                       *values)
+                       Number::size()>              &constraint_mask,
+      Number                                        *values)
     {
       const unsigned int given_degree =
         fe_degree != -1 ? fe_degree : shape_info.data.front().fe_degree;
 
-      const Number *DEAL_II_RESTRICT weights =
+      const Number2 *DEAL_II_RESTRICT weights =
         shape_info.data.front().subface_interpolation_matrices[0].data();
 
       const unsigned int points = given_degree + 1;
@@ -1524,14 +1527,14 @@ namespace internal
     }
 
   public:
-    template <bool transpose>
+    template <bool transpose, typename Number2>
     static void
     run_internal(
-      const unsigned int                            n_desired_components,
-      const MatrixFreeFunctions::ShapeInfo<Number> &shape_info,
+      const unsigned int                             n_desired_components,
+      const MatrixFreeFunctions::ShapeInfo<Number2> &shape_info,
       const std::array<MatrixFreeFunctions::compressed_constraint_kind,
-                       Number::size()>             &constraint_mask,
-      Number                                       *values)
+                       Number::size()>              &constraint_mask,
+      Number                                        *values)
     {
       const unsigned int given_degree =
         fe_degree != -1 ? fe_degree : shape_info.data.front().fe_degree;
@@ -1703,14 +1706,14 @@ namespace internal
   struct FEEvaluationImplHangingNodes
   {
   public:
-    template <int fe_degree>
+    template <int fe_degree, typename Number2>
     static bool
-    run(const unsigned int                            n_desired_components,
-        const MatrixFreeFunctions::ShapeInfo<Number> &shape_info,
-        const bool                                    transpose,
+    run(const unsigned int                             n_desired_components,
+        const MatrixFreeFunctions::ShapeInfo<Number2> &shape_info,
+        const bool                                     transpose,
         const std::array<MatrixFreeFunctions::compressed_constraint_kind,
-                         Number::size()>             &c_mask,
-        Number                                       *values)
+                         Number::size()>              &c_mask,
+        Number                                        *values)
     {
       using RunnerType =
         FEEvaluationImplHangingNodesRunner<used_runner_type<fe_degree>(),
index 7f8f7c360342a51e14318127f507311cac968277..b50575eceaf7f2061cfdfd553ece0a1d855c13e9 100644 (file)
@@ -153,13 +153,13 @@ namespace internal
   struct FEEvaluationHangingNodesFactory
   {
     static void
-    apply(const unsigned int n_components,
-          const unsigned int fe_degree,
-          const MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> &shape_info,
-          const bool                                                 transpose,
+    apply(const unsigned int                             n_components,
+          const unsigned int                             fe_degree,
+          const MatrixFreeFunctions::ShapeInfo<Number>  &shape_info,
+          const bool                                     transpose,
           const std::array<MatrixFreeFunctions::compressed_constraint_kind,
-                           VectorizedArrayType::size()>             &c_mask,
-          VectorizedArrayType                                       *values);
+                           VectorizedArrayType::size()> &c_mask,
+          VectorizedArrayType                           *values);
   };
 
 } // end of namespace internal
index bcfcb992bac18bb009676ea93258d2b2d275275b..ed41d8fab460d2b0ebf4644c4990af4af9276f39 100644 (file)
@@ -31,13 +31,13 @@ namespace internal
   template <int dim, typename Number, typename VectorizedArrayType>
   void
   FEEvaluationHangingNodesFactory<dim, Number, VectorizedArrayType>::apply(
-    const unsigned int                                         n_components,
-    const unsigned int                                         fe_degree,
-    const MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> &shape_info,
-    const bool                                                 transpose,
+    const unsigned int                             n_components,
+    const unsigned int                             fe_degree,
+    const MatrixFreeFunctions::ShapeInfo<Number>  &shape_info,
+    const bool                                     transpose,
     const std::array<MatrixFreeFunctions::compressed_constraint_kind,
-                     VectorizedArrayType::size()>             &c_mask,
-    VectorizedArrayType                                       *values)
+                     VectorizedArrayType::size()> &c_mask,
+    VectorizedArrayType                           *values)
   {
     instantiation_helper_degree_run<
       1,
index 58fc2bae7f15c3e82255566b64a279b4db701b56..4db72c6a376bd90455e99d873cc22644393e1516 100644 (file)
@@ -2439,11 +2439,10 @@ inline FEEvaluationBase<dim,
   (void)base_element_number;
 
   Assert(this->data == nullptr, ExcInternalError());
-  this->data =
-    new internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>(
-      Quadrature<(is_face ? dim - 1 : dim)>(quadrature),
-      fe,
-      fe.component_to_base_index(first_selected_component).first);
+  this->data = new internal::MatrixFreeFunctions::ShapeInfo<Number>(
+    Quadrature<(is_face ? dim - 1 : dim)>(quadrature),
+    fe,
+    fe.component_to_base_index(first_selected_component).first);
 
   this->set_data_pointers(scratch_data_array, n_components_);
 }
@@ -2475,8 +2474,7 @@ inline FEEvaluationBase<dim,
     {
       Assert(other.mapped_geometry.get() != nullptr, ExcInternalError());
       this->data =
-        new internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>(
-          *other.data);
+        new internal::MatrixFreeFunctions::ShapeInfo<Number>(*other.data);
 
       // Create deep copy of mapped geometry for use in parallel
       this->mapped_geometry =
@@ -2543,8 +2541,7 @@ operator=(const FEEvaluationBase<dim,
     {
       Assert(other.mapped_geometry.get() != nullptr, ExcInternalError());
       this->data =
-        new internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>(
-          *other.data);
+        new internal::MatrixFreeFunctions::ShapeInfo<Number>(*other.data);
       scratch_data_array = new AlignedVector<VectorizedArrayType>();
 
       // Create deep copy of mapped geometry for use in parallel
index f9dcc3b8d7dac1d40705e50f4e23120ad526e6aa..9e955ee19b6f67051c759c2ebcc103d95abd53fe 100644 (file)
@@ -112,7 +112,8 @@ namespace internal
 template <int dim, typename Number, bool is_face>
 class FEEvaluationData
 {
-  using ShapeInfoType = internal::MatrixFreeFunctions::ShapeInfo<Number>;
+  using ShapeInfoType = internal::MatrixFreeFunctions::ShapeInfo<
+    typename internal::VectorizedArrayTrait<Number>::value_type>;
   using MappingInfoStorageType = internal::MatrixFreeFunctions::
     MappingInfoStorage<(is_face ? dim - 1 : dim), dim, Number>;
   using DoFInfo = internal::MatrixFreeFunctions::DoFInfo;
@@ -120,10 +121,10 @@ class FEEvaluationData
 public:
   static constexpr unsigned int dimension = dim;
 
-  using NumberType             = Number;
-  using shape_info_number_type = Number;
+  using NumberType = Number;
   using ScalarNumber =
     typename internal::VectorizedArrayTrait<Number>::value_type;
+  using shape_info_number_type = ScalarNumber;
 
   static constexpr unsigned int n_lanes =
     internal::VectorizedArrayTrait<Number>::width();
index 2d11f50e2f6c2681fbc23dfbc6f30211735b3822..44d38e9c6e06a793e9cf50d0e9b9cb5a28d3eba5 100644 (file)
@@ -1151,13 +1151,13 @@ namespace internal
                 typename VectorizedDouble>
       void
       mapping_q_compute_range(
-        const unsigned int                 begin_cell,
-        const unsigned int                 end_cell,
-        const std::vector<GeometryType>   &cell_type,
-        const std::vector<bool>           &process_cell,
-        const UpdateFlags                  update_flags_cells,
-        const AlignedVector<double>       &plain_quadrature_points,
-        const ShapeInfo<VectorizedDouble> &shape_info,
+        const unsigned int               begin_cell,
+        const unsigned int               end_cell,
+        const std::vector<GeometryType> &cell_type,
+        const std::vector<bool>         &process_cell,
+        const UpdateFlags                update_flags_cells,
+        const AlignedVector<double>     &plain_quadrature_points,
+        const ShapeInfo<double>         &shape_info,
         MappingInfoStorage<dim, dim, VectorizedArrayType> &my_data)
       {
         constexpr unsigned int n_lanes   = VectorizedArrayType::size();
@@ -2083,12 +2083,12 @@ namespace internal
         const unsigned int begin_face,
         const unsigned int end_face,
         const std::vector<FaceToCellTopology<VectorizedArrayType::size()>>
-                                          &faces,
-        const std::vector<GeometryType>   &face_type,
-        const std::vector<bool>           &process_face,
-        const UpdateFlags                  update_flags_faces,
-        const AlignedVector<double>       &plain_quadrature_points,
-        const ShapeInfo<VectorizedDouble> &shape_info,
+                                        &faces,
+        const std::vector<GeometryType> &face_type,
+        const std::vector<bool>         &process_face,
+        const UpdateFlags                update_flags_faces,
+        const AlignedVector<double>     &plain_quadrature_points,
+        const ShapeInfo<double>         &shape_info,
         MappingInfoStorage<dim - 1, dim, VectorizedArrayType> &my_data)
       {
         constexpr unsigned int n_lanes   = VectorizedArrayType::size();
@@ -2696,7 +2696,7 @@ namespace internal
       // functions or the quadrature points; shape info is merely a vehicle to
       // return us the right interpolation matrices from the cell support
       // points to the cell and face quadrature points.
-      std::vector<ShapeInfo<VectorizedDouble>> shape_infos(cell_data.size());
+      std::vector<ShapeInfo<double>> shape_infos(cell_data.size());
       {
         FE_DGQ<dim> fe_geometry(mapping_degree);
         for (unsigned int my_q = 0; my_q < cell_data.size(); ++my_q)
index df1deeb17d8338dd10fcb272699f84f970b73525..78410ba48dd2cdb2e0f37e66aa4ea88d85cd2eb6 100644 (file)
@@ -2147,7 +2147,7 @@ public:
   /**
    * Return the unit cell information for given hp-index.
    */
-  const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> &
+  const internal::MatrixFreeFunctions::ShapeInfo<Number> &
   get_shape_info(const unsigned int dof_handler_index_component = 0,
                  const unsigned int quad_index                  = 0,
                  const unsigned int fe_base_element             = 0,
@@ -2292,8 +2292,7 @@ private:
   /**
    * Contains shape value information on the unit cell.
    */
-  Table<4, internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>>
-    shape_info;
+  Table<4, internal::MatrixFreeFunctions::ShapeInfo<Number>> shape_info;
 
   /**
    * Describes how the cells are gone through. With the cell level (first
@@ -2856,7 +2855,7 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_ghost_set(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-inline const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> &
+inline const internal::MatrixFreeFunctions::ShapeInfo<Number> &
 MatrixFree<dim, Number, VectorizedArrayType>::get_shape_info(
   const unsigned int dof_handler_index,
   const unsigned int index_quad,
@@ -3920,7 +3919,17 @@ namespace internal
     zero_vector_region(const unsigned int range_index, VectorType &vec) const
     {
       if (range_index == numbers::invalid_unsigned_int || range_index == 0)
-        vec = typename VectorType::value_type();
+        {
+          if constexpr (std::is_same_v<
+                          ArrayView<typename VectorType::value_type>,
+                          VectorType>)
+            {
+              for (unsigned int i = 0; i < vec.size(); ++i)
+                vec[i] = typename VectorType::value_type();
+            }
+          else
+            vec = typename VectorType::value_type();
+        }
     }
 
 
index abaddb0fdc3f8feba84e622617da8f9c2abb3d36..ae41d860663d006dde8fe847c7ef79c23df52419 100644 (file)
@@ -74,19 +74,18 @@ namespace dealii
 
         };
 
-    const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
-      *shape_info_base;
+    const internal::MatrixFreeFunctions::ShapeInfo<Number> *shape_info_base;
 
     void
     reinit(unsigned int cell_batch_index,
-           const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
-             *shape_info = nullptr);
+           const internal::MatrixFreeFunctions::ShapeInfo<Number> *shape_info =
+             nullptr);
 
     void
     reinit(unsigned int cell_batch_index,
            unsigned int face_number,
-           const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
-             *shape_info = nullptr);
+           const internal::MatrixFreeFunctions::ShapeInfo<Number> *shape_info =
+             nullptr);
   };
 
   template <int dim,
@@ -103,8 +102,7 @@ namespace dealii
                     Number,
                     VectorizedArrayType>::
     reinit(const unsigned int cell_batch_index,
-           const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
-             *shape_info)
+           const internal::MatrixFreeFunctions::ShapeInfo<Number> *shape_info)
   {
     Assert(this->mapped_geometry == nullptr,
            ExcMessage(
@@ -163,8 +161,7 @@ namespace dealii
                     VectorizedArrayType>::
     reinit(const unsigned int cell_batch_index,
            const unsigned int face_number,
-           const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
-             *shape_info)
+           const internal::MatrixFreeFunctions::ShapeInfo<Number> *shape_info)
   {
     Assert(
       this->quad_no <
@@ -344,8 +341,7 @@ test(const dealii::FE_Poly<dim> &fe)
 
   matrix_free.initialize_dof_vector(dst);
 
-  std::array<internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>,
-             2 * dim>
+  std::array<internal::MatrixFreeFunctions::ShapeInfo<Number>, 2 * dim>
     shape_info_shift;
 
   const QGauss<1> quadrature_1D(n_points);
@@ -355,11 +351,11 @@ test(const dealii::FE_Poly<dim> &fe)
   const Quadrature<1> quadrature_1D_shift_plus =
     shift_1d_quadrature(quadrature_1D, +1.);
 
-  dealii::internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
-    shape_info_base(quadrature_1D, fe);
-  dealii::internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
+  dealii::internal::MatrixFreeFunctions::ShapeInfo<Number> shape_info_base(
+    quadrature_1D, fe);
+  dealii::internal::MatrixFreeFunctions::ShapeInfo<Number>
     shape_info_shift_plus(quadrature_1D_shift_plus, fe);
-  dealii::internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType>
+  dealii::internal::MatrixFreeFunctions::ShapeInfo<Number>
     shape_info_shift_minus(quadrature_1D_shift_minus, fe);
 
   for (unsigned int f = 0; f < 2 * dim; ++f)
index 49a6520476873895bd0e87eb35cd9f8e3f36ffca..fc6673718713d55a12ca05f418cf0843863e70bc 100644 (file)
@@ -202,10 +202,10 @@ namespace dealii
                               const auto w =
                                 transpose ?
                                   constraint_weights[i * (fe_degree + 1) +
-                                                     interp_idx][v] :
+                                                     interp_idx] :
                                   constraint_weights[interp_idx *
                                                        (fe_degree + 1) +
-                                                     i][v];
+                                                     i];
                               t += w * values_temp[real_idx][v];
                             }
                         }
@@ -222,11 +222,10 @@ namespace dealii
                                 transpose ?
                                   constraint_weights[(fe_degree - i) *
                                                        (fe_degree + 1) +
-                                                     fe_degree - interp_idx]
-                                                    [v] :
+                                                     fe_degree - interp_idx] :
                                   constraint_weights[(fe_degree - interp_idx) *
                                                        (fe_degree + 1) +
-                                                     fe_degree - i][v];
+                                                     fe_degree - i];
                               t += w * values_temp[real_idx][v];
                             }
                         }
@@ -375,10 +374,10 @@ namespace dealii
                                 const auto w =
                                   transpose ?
                                     constraint_weights[i * (fe_degree + 1) +
-                                                       interp_idx][v] :
+                                                       interp_idx] :
                                     constraint_weights[interp_idx *
                                                          (fe_degree + 1) +
-                                                       i][v];
+                                                       i];
                                 t += w * values_temp[real_idx][v];
                               }
                           }
@@ -397,12 +396,11 @@ namespace dealii
                                   transpose ?
                                     constraint_weights[(fe_degree - i) *
                                                          (fe_degree + 1) +
-                                                       fe_degree - interp_idx]
-                                                      [v] :
+                                                       fe_degree - interp_idx] :
                                     constraint_weights[(fe_degree -
                                                         interp_idx) *
                                                          (fe_degree + 1) +
-                                                       fe_degree - i][v];
+                                                       fe_degree - i];
                                 t += w * values_temp[real_idx][v];
                               }
                           }
index bf322a51f99bdc7902696262bf67e29e6b1ef00f..849779097bdcfd7b576a8da20d62b049305de70a 100644 (file)
@@ -674,17 +674,19 @@ namespace NavierStokes_DG
       dim,
       n_points_1d,
       n_points_1d,
-      VectorizedArrayType>
-      eval(AlignedVector<VectorizedArrayType>(),
+      VectorizedArrayType,
+      Number>
+      eval({},
            data.get_shape_info().data[0].shape_gradients_collocation_eo,
-           AlignedVector<VectorizedArrayType>());
+           {});
 
     internal::EvaluatorTensorProduct<
       internal::EvaluatorVariant::evaluate_evenodd,
       dim - 1,
       n_q_points_1d,
       n_q_points_1d,
-      VectorizedArrayType>
+      VectorizedArrayType,
+      Number>
       eval_face({},
                 data.get_shape_info().data[0].shape_gradients_collocation_eo,
                 {});
@@ -789,7 +791,7 @@ namespace NavierStokes_DG
 
             phi_m.reinit(cell, face);
 
-            const AlignedVector<VectorizedArrayType> &shape_data =
+            const AlignedVector<Number> &shape_data =
               data.get_shape_info().data.front().quadrature_data_on_face[face %
                                                                          2];
             const std::array<int, 2> n_blocks{

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.