]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Shorten compilation time of two MatrixFree tests
authorMartin Kronbichler <martin.kronbichler@uni-a.de>
Thu, 29 Jun 2023 10:34:31 +0000 (12:34 +0200)
committerMartin Kronbichler <martin.kronbichler@uni-a.de>
Fri, 30 Jun 2023 05:57:36 +0000 (07:57 +0200)
tests/matrix_free/matrix_vector_hp.cc
tests/matrix_free/matrix_vector_hp.output
tests/matrix_free/matrix_vector_mf.h
tests/matrix_free/thread_correctness_hp.cc

index fbba9210b572bc19d9fb10ff5f5905cdfb915008..8dcc984c3912c1084e43398b87ffe265ac075ee0 100644 (file)
 #include <deal.II/base/function.h>
 
 #include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/manifold_lib.h>
 
 #include <deal.II/hp/fe_values.h>
 
+#include <deal.II/lac/dynamic_sparsity_pattern.h>
+#include <deal.II/lac/sparse_matrix.h>
+
+#include <deal.II/numerics/vector_tools.h>
+
 #include "../tests.h"
 
-#include "matrix_vector_common.h"
+#include "matrix_vector_mf.h"
 
 
 
@@ -45,7 +54,6 @@ public:
               const Vector<Number> &                       src,
               const std::pair<unsigned int, unsigned int> &cell_range) const
   {
-    // ask MatrixFree for cell_range for different orders
     std::pair<unsigned int, unsigned int> subrange_deg =
       data.create_cell_subrange_hp(cell_range, 1);
     if (subrange_deg.second > subrange_deg.first)
@@ -65,30 +73,6 @@ public:
                                                     dst,
                                                     src,
                                                     subrange_deg);
-    subrange_deg = data.create_cell_subrange_hp(cell_range, 4);
-    if (subrange_deg.second > subrange_deg.first)
-      helmholtz_operator<dim, 4, Vector<Number>, 5>(data,
-                                                    dst,
-                                                    src,
-                                                    subrange_deg);
-    subrange_deg = data.create_cell_subrange_hp(cell_range, 5);
-    if (subrange_deg.second > subrange_deg.first)
-      helmholtz_operator<dim, 5, Vector<Number>, 6>(data,
-                                                    dst,
-                                                    src,
-                                                    subrange_deg);
-    subrange_deg = data.create_cell_subrange_hp(cell_range, 6);
-    if (subrange_deg.second > subrange_deg.first)
-      helmholtz_operator<dim, 6, Vector<Number>, 7>(data,
-                                                    dst,
-                                                    src,
-                                                    subrange_deg);
-    subrange_deg = data.create_cell_subrange_hp(cell_range, 7);
-    if (subrange_deg.second > subrange_deg.first)
-      helmholtz_operator<dim, 7, Vector<Number>, 8>(data,
-                                                    dst,
-                                                    src,
-                                                    subrange_deg);
   }
 
   void
@@ -104,20 +88,15 @@ private:
 
 
 
-template <int dim, int fe_degree>
+template <int dim>
 void
 test()
 {
-  if (fe_degree > 1)
-    return;
-
   using number = double;
   const SphericalManifold<dim> manifold;
   Triangulation<dim>           tria;
   GridGenerator::hyper_ball(tria);
-  typename Triangulation<dim>::active_cell_iterator cell = tria.begin_active(),
-                                                    endc = tria.end();
-  for (; cell != endc; ++cell)
+  for (const auto &cell : tria.active_cell_iterators())
     for (const unsigned int f : GeometryInfo<dim>::face_indices())
       if (cell->at_boundary(f))
         cell->face(f)->set_all_manifold_ids(0);
@@ -127,17 +106,14 @@ test()
   // refine a few cells
   for (unsigned int i = 0; i < 11 - 3 * dim; ++i)
     {
-      typename Triangulation<dim>::active_cell_iterator cell =
-                                                          tria.begin_active(),
-                                                        endc = tria.end();
-      unsigned int counter                                   = 0;
-      for (; cell != endc; ++cell, ++counter)
-        if (counter % (7 - i) == 0)
+      unsigned int counter = 0;
+      for (const auto &cell : tria.active_cell_iterators())
+        if ((counter++) % (7 - i) == 0)
           cell->set_refine_flag();
       tria.execute_coarsening_and_refinement();
     }
 
-  const unsigned int max_degree = 9 - 2 * dim;
+  const unsigned int max_degree = 3;
 
   hp::FECollection<dim> fe_collection;
   hp::QCollection<dim>  quadrature_collection;
@@ -156,9 +132,7 @@ test()
   DoFHandler<dim> dof(tria);
   // set the active FE index in a random order
   {
-    typename DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
-                                                   endc = dof.end();
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof.active_cell_iterators())
       {
         const unsigned int fe_index = Testing::rand() % max_degree;
         cell->set_active_fe_index(fe_index + 1);
@@ -166,8 +140,7 @@ test()
     // We cannot set random cells to FE_Nothing. We get the following error
     // The violated condition was: dominating_fe.n_dofs_per_face(face) <=
     // subface_fe.n_dofs_per_face(face)
-    cell = dof.begin_active();
-    cell->set_active_fe_index(0);
+    dof.begin_active()->set_active_fe_index(0);
   }
 
   // setup DoFs
@@ -208,9 +181,7 @@ test()
     FullMatrix<double>                   cell_matrix;
     std::vector<types::global_dof_index> local_dof_indices;
 
-    typename DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
-                                                   endc = dof.end();
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof.active_cell_iterators())
       {
         const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell;
 
@@ -259,3 +230,13 @@ test()
   const double diff_norm = result_mf.linfty_norm();
   deallog << "Norm of difference: " << diff_norm << std::endl << std::endl;
 }
+
+
+
+int
+main()
+{
+  initlog();
+
+  test<2>();
+}
index 7252835322c1b892592e6910f887cad1cad53fe3..444fbc53b0eed5d39b445cda156266c58b2d171e 100644 (file)
@@ -1,5 +1,3 @@
 
-DEAL:2d::Norm of difference: 0
-DEAL:2d::
-DEAL:3d::Norm of difference: 0
-DEAL:3d::
+DEAL::Norm of difference: 1.06581e-14
+DEAL::
index 984515067a88428666dc3da873a0e9ed4455e9f5..877dd5efd6888107806db7b1568ea3c95f634aab 100644 (file)
@@ -56,6 +56,35 @@ helmholtz_operator(const MatrixFree<dim, typename VectorType::value_type> &data,
 
 
 
+template <int dim, typename VectorType>
+void
+helmholtz_operator_no_template(
+  const MatrixFree<dim, typename VectorType::value_type> &data,
+  VectorType &                                            dst,
+  const VectorType &                                      src,
+  const std::pair<unsigned int, unsigned int> &           cell_range)
+{
+  using Number = typename VectorType::value_type;
+  FEEvaluation<dim, -1, 0, 1, Number> fe_eval(data, cell_range);
+  const unsigned int                  n_q_points = fe_eval.n_q_points;
+
+  for (unsigned int cell = cell_range.first; cell < cell_range.second; ++cell)
+    {
+      fe_eval.reinit(cell);
+      fe_eval.read_dof_values(src);
+      fe_eval.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
+      for (unsigned int q = 0; q < n_q_points; ++q)
+        {
+          fe_eval.submit_value(Number(10) * fe_eval.get_value(q), q);
+          fe_eval.submit_gradient(fe_eval.get_gradient(q), q);
+        }
+      fe_eval.integrate(EvaluationFlags::values | EvaluationFlags::gradients);
+      fe_eval.distribute_local_to_global(dst);
+    }
+}
+
+
+
 template <int dim, typename VectorType>
 void
 helmholtz_operator_no_template(
index 88a6d7df4cf3d46bd28d48a9f1c28b407e719bbf..c85f45810080fe04598fa2b56bb7922700f7c568 100644 (file)
@@ -46,26 +46,7 @@ public:
               const Vector<Number> &                       src,
               const std::pair<unsigned int, unsigned int> &cell_range) const
   {
-    // Ask MatrixFree for cell_range for different
-    // orders
-    std::pair<unsigned int, unsigned int> subrange_deg;
-#define CALL_METHOD(degree)                                         \
-  subrange_deg = data.create_cell_subrange_hp(cell_range, degree);  \
-  if (subrange_deg.second > subrange_deg.first)                     \
-  helmholtz_operator<dim, degree, Vector<Number>, degree + 1>(data, \
-                                                              dst,  \
-                                                              src,  \
-                                                              subrange_deg)
-
-    CALL_METHOD(1);
-    CALL_METHOD(2);
-    CALL_METHOD(3);
-    CALL_METHOD(4);
-    CALL_METHOD(5);
-    CALL_METHOD(6);
-    CALL_METHOD(7);
-
-#undef CALL_METHOD
+    helmholtz_operator_no_template<dim>(data, dst, src, cell_range);
   }
 
   void

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.