]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clean up weight_fe_q_dofs_by_entity 14302/head
authorPeter Munch <peterrmuench@gmail.com>
Thu, 22 Sep 2022 08:24:40 +0000 (10:24 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Thu, 22 Sep 2022 11:06:36 +0000 (13:06 +0200)
include/deal.II/matrix_free/tensor_product_kernels.h
include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h
source/multigrid/mg_transfer_matrix_free.cc

index 88b0a8b2d5197d2377d603f959c4e032dd44c8a2..32cd53edb2605473ec97457e550462d314ed5d5b 100644 (file)
@@ -3324,36 +3324,39 @@ namespace internal
   }
 
 
-  template <int dim, int loop_length_template, typename Number>
+  template <int dim, int n_points_1d_template, typename Number>
   inline void
-  weight_fe_q_dofs_by_entity(const VectorizedArray<Number> *weights,
-                             const unsigned int             n_components,
-                             const int                loop_length_non_template,
-                             VectorizedArray<Number> *data)
+  weight_fe_q_dofs_by_entity(const Number *     weights,
+                             const unsigned int n_components,
+                             const int          n_points_1d_non_template,
+                             Number *           data)
   {
-    const int loop_length = loop_length_template != -1 ?
-                              loop_length_template :
-                              loop_length_non_template;
-
-    Assert(loop_length > 0, ExcNotImplemented());
-    Assert(loop_length < 100, ExcNotImplemented());
-    unsigned int degree_to_3[100];
-    degree_to_3[0] = 0;
-    for (int i = 1; i < loop_length - 1; ++i)
-      degree_to_3[i] = 1;
-    degree_to_3[loop_length - 1] = 2;
+    const int n_points_1d = n_points_1d_template != -1 ?
+                              n_points_1d_template :
+                              n_points_1d_non_template;
+
+    Assert(n_points_1d > 0, ExcNotImplemented());
+    Assert(n_points_1d < 100, ExcNotImplemented());
+
+    unsigned int compressed_index[100];
+    compressed_index[0] = 0;
+    for (int i = 1; i < n_points_1d - 1; ++i)
+      compressed_index[i] = 1;
+    compressed_index[n_points_1d - 1] = 2;
+
     for (unsigned int c = 0; c < n_components; ++c)
-      for (int k = 0; k < (dim > 2 ? loop_length : 1); ++k)
-        for (int j = 0; j < (dim > 1 ? loop_length : 1); ++j)
+      for (int k = 0; k < (dim > 2 ? n_points_1d : 1); ++k)
+        for (int j = 0; j < (dim > 1 ? n_points_1d : 1); ++j)
           {
-            const unsigned int shift = 9 * degree_to_3[k] + 3 * degree_to_3[j];
+            const unsigned int shift =
+              9 * compressed_index[k] + 3 * compressed_index[j];
             data[0] *= weights[shift];
-            // loop bound as int avoids compiler warnings in case loop_length
+            // loop bound as int avoids compiler warnings in case n_points_1d
             // == 1 (polynomial degree 0)
-            for (int i = 1; i < loop_length - 1; ++i)
+            for (int i = 1; i < n_points_1d - 1; ++i)
               data[i] *= weights[shift + 1];
-            data[loop_length - 1] *= weights[shift + 2];
-            data += loop_length;
+            data[n_points_1d - 1] *= weights[shift + 2];
+            data += n_points_1d;
           }
   }
 
index 08de0842e91ce660f4914efcc14ed0cdb3a152f1..ce52dda65c15d65b5bf973e59e593e1bdfdf80e0 100644 (file)
@@ -2590,11 +2590,12 @@ MGTwoLevelTransfer<dim, LinearAlgebra::distributed::Vector<Number>>::
           // weight and write into dst vector
           if (fine_element_is_continuous && this->weights_compressed.size() > 0)
             {
-              internal::weight_fe_q_dofs_by_entity<dim, -1, Number>(
-                weights_compressed,
-                n_components,
-                scheme.degree_fine + 1,
-                evaluation_data_fine.begin());
+              internal::
+                weight_fe_q_dofs_by_entity<dim, -1, VectorizedArray<Number>>(
+                  weights_compressed,
+                  n_components,
+                  scheme.degree_fine + 1,
+                  evaluation_data_fine.begin());
               weights_compressed += Utilities::pow(3, dim);
             }
 
@@ -2725,11 +2726,12 @@ MGTwoLevelTransfer<dim, LinearAlgebra::distributed::Vector<Number>>::
 
           if (fine_element_is_continuous && this->weights_compressed.size() > 0)
             {
-              internal::weight_fe_q_dofs_by_entity<dim, -1, Number>(
-                weights_compressed,
-                n_components,
-                scheme.degree_fine + 1,
-                evaluation_data_fine.begin());
+              internal::
+                weight_fe_q_dofs_by_entity<dim, -1, VectorizedArray<Number>>(
+                  weights_compressed,
+                  n_components,
+                  scheme.degree_fine + 1,
+                  evaluation_data_fine.begin());
               weights_compressed += Utilities::pow(3, dim);
             }
 
index 6db4b3638bf750e1c051463e649a0b2141ca9381..a131d8846b39b2b5bee81c5055448411e8ff54ec 100644 (file)
@@ -511,7 +511,7 @@ MGTransferMatrixFree<dim, Number>::do_prolongate_add(
           internal::weight_fe_q_dofs_by_entity<dim,
                                                degree != -1 ? 2 * degree + 1 :
                                                               -1,
-                                               Number>(
+                                               VectorizedArray<Number>>(
             &weights_on_refined[to_level - 1][(cell / vec_size) * three_to_dim],
             n_components,
             2 * fe_degree + 1,
@@ -592,14 +592,15 @@ MGTransferMatrixFree<dim, Number>::do_restrict_add(
       // perform tensorized operation
       if (element_is_continuous)
         {
-          internal::weight_fe_q_dofs_by_entity<
-            dim,
-            degree != -1 ? 2 * degree + 1 : -1,
-            Number>(&weights_on_refined[from_level - 1]
-                                       [(cell / vec_size) * three_to_dim],
-                    n_components,
-                    2 * fe_degree + 1,
-                    evaluation_data.data());
+          internal::weight_fe_q_dofs_by_entity<dim,
+                                               degree != -1 ? 2 * degree + 1 :
+                                                              -1,
+                                               VectorizedArray<Number>>(
+            &weights_on_refined[from_level - 1]
+                               [(cell / vec_size) * three_to_dim],
+            n_components,
+            2 * fe_degree + 1,
+            evaluation_data.data());
           for (unsigned int c = 0; c < n_components; ++c)
             internal::FEEvaluationImplBasisChange<
               internal::evaluate_general,

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.