]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Split up matrix-free compilation into more files 10003/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 1 May 2020 07:38:39 +0000 (09:38 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 1 May 2020 09:34:19 +0000 (11:34 +0200)
13 files changed:
include/deal.II/matrix_free/dof_info.h
include/deal.II/matrix_free/dof_info.templates.h
include/deal.II/matrix_free/face_setup_internal.h
include/deal.II/matrix_free/matrix_free.templates.h
include/deal.II/matrix_free/task_info.h
source/matrix_free/CMakeLists.txt
source/matrix_free/dof_info.cc [new file with mode: 0644]
source/matrix_free/matrix_free.cc
source/matrix_free/matrix_free.inst.in
source/matrix_free/matrix_free_inst2.cc [new file with mode: 0644]
source/matrix_free/matrix_free_inst3.cc [new file with mode: 0644]
source/matrix_free/shape_info.cc [new file with mode: 0644]
source/matrix_free/shape_info.inst.in [new file with mode: 0644]

index 7da0f1bf1c68a5732d670ca5f599d04e44300d25..96d6286d1f85f88801933db0f69465936fd2eb84 100644 (file)
@@ -30,6 +30,7 @@
 #include <deal.II/lac/dynamic_sparsity_pattern.h>
 
 #include <deal.II/matrix_free/face_info.h>
+#include <deal.II/matrix_free/mapping_info.h>
 #include <deal.II/matrix_free/task_info.h>
 
 #include <array>
@@ -42,11 +43,38 @@ namespace internal
 {
   namespace MatrixFreeFunctions
   {
-    // Forward declaration
-#ifndef DOXYGEN
+    /**
+     * A struct that takes entries describing a constraint and puts them into
+     * a sorted list where duplicates are filtered out
+     */
     template <typename Number>
-    struct ConstraintValues;
-#endif
+    struct ConstraintValues
+    {
+      ConstraintValues();
+
+      /**
+       * This function inserts some constrained entries to the collection of
+       * all values. It stores the (reordered) numbering of the dofs
+       * (according to the ordering that matches with the function) in
+       * new_indices, and returns the storage position the double array for
+       * access later on.
+       */
+      template <typename number2>
+      unsigned short
+      insert_entries(
+        const std::vector<std::pair<types::global_dof_index, number2>>
+          &entries);
+
+      std::vector<std::pair<types::global_dof_index, double>>
+                                           constraint_entries;
+      std::vector<types::global_dof_index> constraint_indices;
+
+      std::pair<std::vector<Number>, types::global_dof_index> next_constraint;
+      std::map<std::vector<Number>,
+               types::global_dof_index,
+               FPArrayComparator<Number>>
+        constraints;
+    };
 
     /**
      * The class that stores the indices of the degrees of freedom for all the
index 724197e3c6b21f1e5e5bb2144f1b035b3e4c1e86..e2dd2915107955f95a91b3a2b1a9b36625298284 100644 (file)
@@ -37,55 +37,13 @@ namespace internal
 {
   namespace MatrixFreeFunctions
   {
-    struct ConstraintComparator
-    {
-      bool
-      operator()(const std::pair<types::global_dof_index, double> &p1,
-                 const std::pair<types::global_dof_index, double> &p2) const
-      {
-        return p1.second < p2.second;
-      }
-    };
-
-    /**
-     * A struct that takes entries describing a constraint and puts them into
-     * a sorted list where duplicates are filtered out
-     */
-    template <typename Number>
-    struct ConstraintValues
-    {
-      ConstraintValues();
-
-      /**
-       * This function inserts some constrained entries to the collection of
-       * all values. It stores the (reordered) numbering of the dofs
-       * (according to the ordering that matches with the function) in
-       * new_indices, and returns the storage position the double array for
-       * access later on.
-       */
-      template <typename number2>
-      unsigned short
-      insert_entries(
-        const std::vector<std::pair<types::global_dof_index, number2>>
-          &entries);
-
-      std::vector<std::pair<types::global_dof_index, double>>
-                                           constraint_entries;
-      std::vector<types::global_dof_index> constraint_indices;
-
-      std::pair<std::vector<Number>, types::global_dof_index> next_constraint;
-      std::map<std::vector<Number>,
-               types::global_dof_index,
-               FPArrayComparator<double>>
-        constraints;
-    };
-
-
     template <typename Number>
     ConstraintValues<Number>::ConstraintValues()
-      : constraints(FPArrayComparator<double>(1.))
+      : constraints(FPArrayComparator<Number>(1.))
     {}
 
+
+
     template <typename Number>
     template <typename number2>
     unsigned short
@@ -101,7 +59,10 @@ namespace internal
           constraint_entries.assign(entries.begin(), entries.end());
           std::sort(constraint_entries.begin(),
                     constraint_entries.end(),
-                    ConstraintComparator());
+                    [](const std::pair<types::global_dof_index, double> &p1,
+                       const std::pair<types::global_dof_index, double> &p2) {
+                      return p1.second < p2.second;
+                    });
           for (types::global_dof_index j = 0; j < constraint_entries.size();
                j++)
             {
@@ -121,11 +82,7 @@ namespace internal
       // equal. this was quite lengthy and now we use a std::map with a
       // user-defined comparator to compare floating point arrays to a
       // tolerance 1e-13.
-      std::pair<typename std::map<std::vector<double>,
-                                  types::global_dof_index,
-                                  FPArrayComparator<double>>::iterator,
-                bool>
-        it = constraints.insert(next_constraint);
+      const auto it = constraints.insert(next_constraint);
 
       types::global_dof_index insert_position = numbers::invalid_dof_index;
       if (it.second == false)
index b019cd3cc4847a0bbf41be0663db33174891aacf..c6cf4465d1306bdd862baf5d39e4c770b0e1406d 100644 (file)
@@ -1010,7 +1010,7 @@ namespace internal
      * face number, subface index and orientation are the same. This is used
      * to batch similar faces together for vectorization.
      */
-    bool
+    inline bool
     compare_faces_for_vectorization(const FaceToCellTopology<1> &face1,
                                     const FaceToCellTopology<1> &face2)
     {
index 783378cdeeb27a198a472cb517d18bb1a4a84cc9..6912b215b9cabb790ff50aa8a750547950a180eb 100644 (file)
 
 #include <deal.II/dofs/dof_accessor.h>
 
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_dgq.h>
 #include <deal.II/fe/fe_poly.h>
+#include <deal.II/fe/fe_q_dg0.h>
 
 #include <deal.II/hp/q_collection.h>
 
-#include <deal.II/matrix_free/dof_info.templates.h>
 #include <deal.II/matrix_free/face_info.h>
 #include <deal.II/matrix_free/face_setup_internal.h>
 #include <deal.II/matrix_free/matrix_free.h>
-#include <deal.II/matrix_free/shape_info.templates.h>
 
 #ifdef DEAL_II_WITH_THREADS
 #  include <deal.II/base/parallel.h>
@@ -1248,19 +1249,14 @@ MatrixFree<dim, Number, VectorizedArrayType>::initialize_indices(
   }
 
   // set constraint pool from the std::map and reorder the indices
-  typename std::map<std::vector<double>,
-                    types::global_dof_index,
-                    internal::MatrixFreeFunctions::FPArrayComparator<double>>::
-    iterator it  = constraint_values.constraints.begin(),
-             end = constraint_values.constraints.end();
   std::vector<const std::vector<double> *> constraints(
     constraint_values.constraints.size());
   unsigned int length = 0;
-  for (; it != end; ++it)
+  for (const auto &it : constraint_values.constraints)
     {
-      AssertIndexRange(it->second, constraints.size());
-      constraints[it->second] = &it->first;
-      length += it->first.size();
+      AssertIndexRange(it.second, constraints.size());
+      constraints[it.second] = &it.first;
+      length += it.first.size();
     }
   constraint_pool_data.clear();
   constraint_pool_data.reserve(length);
@@ -1813,7 +1809,7 @@ MatrixFree<dim, Number, VectorizedArrayType>::clear()
 
 namespace internal
 {
-  void
+  inline void
   fill_index_subrange(
     const unsigned int                                        begin,
     const unsigned int                                        end,
@@ -1832,7 +1828,7 @@ namespace internal
   }
 
   template <int dim>
-  void
+  inline void
   fill_connectivity_subrange(
     const unsigned int                                        begin,
     const unsigned int                                        end,
@@ -1875,7 +1871,7 @@ namespace internal
       }
   }
 
-  void
+  inline void
   fill_connectivity_indirect_subrange(
     const unsigned int            begin,
     const unsigned int            end,
index 12ad101c93eb34adbbdc53c163ce28f39ed1ffb8..cfc972ccd99d63338e504308dddad3ffed7340c0 100644 (file)
@@ -95,10 +95,6 @@ namespace internal
 
   namespace MatrixFreeFunctions
   {
-    // forward declaration of internal data structure
-    template <typename Number>
-    struct ConstraintValues;
-
     /**
      * A struct that collects all information related to parallelization with
      * threads: The work is subdivided into tasks that can be done
index fa27962bdcb3663a02219826fa1d3c56dacb15dd..8c1bb8aeac0ea5524eed73baefbab1e6bc71203d 100644 (file)
 INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR})
 
 SET(_src
+  dof_info.cc
   evaluation_selector.cc
   mapping_info.cc
   mapping_info_inst2.cc
   mapping_info_inst3.cc
   matrix_free.cc
+  matrix_free_inst2.cc
+  matrix_free_inst3.cc
+  shape_info.cc
   task_info.cc
   )
 
@@ -28,6 +32,7 @@ SET(_inst
   evaluation_selector.inst.in
   mapping_info.inst.in
   matrix_free.inst.in
+  shape_info.inst.in
   )
 
 FILE(GLOB _header
diff --git a/source/matrix_free/dof_info.cc b/source/matrix_free/dof_info.cc
new file mode 100644 (file)
index 0000000..105b44d
--- /dev/null
@@ -0,0 +1,109 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/conditional_ostream.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/vectorization.h>
+
+#include <deal.II/matrix_free/dof_info.templates.h>
+
+#include <iostream>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace internal
+{
+  namespace MatrixFreeFunctions
+  {
+    template struct ConstraintValues<double>;
+    template struct ConstraintValues<float>;
+
+    template void
+    DoFInfo::read_dof_indices<double>(
+      const std::vector<types::global_dof_index> &,
+      const std::vector<unsigned int> &,
+      const AffineConstraints<double> &,
+      const unsigned int,
+      ConstraintValues<double> &,
+      bool &);
+
+    template void
+    DoFInfo::read_dof_indices<float>(
+      const std::vector<types::global_dof_index> &,
+      const std::vector<unsigned int> &,
+      const AffineConstraints<float> &,
+      const unsigned int,
+      ConstraintValues<double> &,
+      bool &);
+
+    template void
+    DoFInfo::compute_face_index_compression<1>(
+      const std::vector<FaceToCellTopology<1>> &);
+    template void
+    DoFInfo::compute_face_index_compression<2>(
+      const std::vector<FaceToCellTopology<2>> &);
+    template void
+    DoFInfo::compute_face_index_compression<4>(
+      const std::vector<FaceToCellTopology<4>> &);
+    template void
+    DoFInfo::compute_face_index_compression<8>(
+      const std::vector<FaceToCellTopology<8>> &);
+    template void
+    DoFInfo::compute_face_index_compression<16>(
+      const std::vector<FaceToCellTopology<16>> &);
+
+    template void
+    DoFInfo::compute_vector_zero_access_pattern<1>(
+      const TaskInfo &,
+      const std::vector<FaceToCellTopology<1>> &);
+    template void
+    DoFInfo::compute_vector_zero_access_pattern<2>(
+      const TaskInfo &,
+      const std::vector<FaceToCellTopology<2>> &);
+    template void
+    DoFInfo::compute_vector_zero_access_pattern<4>(
+      const TaskInfo &,
+      const std::vector<FaceToCellTopology<4>> &);
+    template void
+    DoFInfo::compute_vector_zero_access_pattern<8>(
+      const TaskInfo &,
+      const std::vector<FaceToCellTopology<8>> &);
+    template void
+    DoFInfo::compute_vector_zero_access_pattern<16>(
+      const TaskInfo &,
+      const std::vector<FaceToCellTopology<16>> &);
+
+    template void
+    DoFInfo::print_memory_consumption<std::ostream>(std::ostream &,
+                                                    const TaskInfo &) const;
+    template void
+    DoFInfo::print_memory_consumption<ConditionalOStream>(
+      ConditionalOStream &,
+      const TaskInfo &) const;
+
+    template void
+    DoFInfo::print<double>(const std::vector<double> &,
+                           const std::vector<unsigned int> &,
+                           std::ostream &) const;
+
+    template void
+    DoFInfo::print<float>(const std::vector<float> &,
+                          const std::vector<unsigned int> &,
+                          std::ostream &) const;
+  } // namespace MatrixFreeFunctions
+} // namespace internal
+
+DEAL_II_NAMESPACE_CLOSE
index e73c54fd3925bb6b2d1370e19ca7aaeabf640bf4..8f0e9a153a4698610f0be9819475baaf13662c34 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
+#define SPLIT_INSTANTIATIONS_COUNT 3
+#ifndef SPLIT_INSTANTIATIONS_INDEX
+#  define SPLIT_INSTANTIATIONS_INDEX 0
+#endif
 #include "matrix_free.inst"
 
+#if SPLIT_INSTANTIATIONS_INDEX == 0
+
 template struct internal::MatrixFreeFunctions::ShapeInfo<double>;
 template struct internal::MatrixFreeFunctions::ShapeInfo<float>;
 
+#endif
+
 DEAL_II_NAMESPACE_CLOSE
index 22762c607e2f2e6f5403af169644cbb16968931e..53c7631a99a817ca09b578e9bfa1786e00f2efaa 100644 (file)
 // ---------------------------------------------------------------------
 
 
-for (deal_II_dimension : DIMENSIONS; deal_II_scalar : REAL_SCALARS)
-  {
-    template void internal::MatrixFreeFunctions::ShapeInfo<deal_II_scalar>::
-      reinit<deal_II_dimension>(
-        const Quadrature<1> &,
-        const FiniteElement<deal_II_dimension, deal_II_dimension> &,
-        const unsigned int);
-  }
-
-
-
 for (deal_II_dimension : DIMENSIONS;
      deal_II_scalar_vectorized : REAL_SCALARS_VECTORIZED)
   {
@@ -78,12 +67,6 @@ for (deal_II_dimension : DIMENSIONS;
                   deal_II_scalar_vectorized>::
         get_dof_handler<hp::DoFHandler<deal_II_dimension>>(const unsigned int)
           const;
-
-    template void internal::MatrixFreeFunctions::
-      ShapeInfo<deal_II_scalar_vectorized>::reinit<deal_II_dimension>(
-        const Quadrature<1> &,
-        const FiniteElement<deal_II_dimension, deal_II_dimension> &,
-        const unsigned int);
   }
 
 
@@ -125,9 +108,3 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS;
         const FiniteElement<deal_II_dimension, deal_II_space_dimension> &);
 #endif
   }
-
-for (deal_II_scalar_vectorized : REAL_SCALARS_VECTORIZED)
-  {
-    template struct internal::MatrixFreeFunctions::ShapeInfo<
-      deal_II_scalar_vectorized>;
-  }
diff --git a/source/matrix_free/matrix_free_inst2.cc b/source/matrix_free/matrix_free_inst2.cc
new file mode 100644 (file)
index 0000000..998346d
--- /dev/null
@@ -0,0 +1,17 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#define SPLIT_INSTANTIATIONS_INDEX 1
+#include "matrix_free.cc"
diff --git a/source/matrix_free/matrix_free_inst3.cc b/source/matrix_free/matrix_free_inst3.cc
new file mode 100644 (file)
index 0000000..8358b94
--- /dev/null
@@ -0,0 +1,17 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#define SPLIT_INSTANTIATIONS_INDEX 2
+#include "matrix_free.cc"
diff --git a/source/matrix_free/shape_info.cc b/source/matrix_free/shape_info.cc
new file mode 100644 (file)
index 0000000..aed24a6
--- /dev/null
@@ -0,0 +1,32 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/conditional_ostream.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/vectorization.h>
+
+#include <deal.II/matrix_free/shape_info.templates.h>
+
+#include <iostream>
+
+DEAL_II_NAMESPACE_OPEN
+
+#include "shape_info.inst"
+
+template struct internal::MatrixFreeFunctions::ShapeInfo<double>;
+template struct internal::MatrixFreeFunctions::ShapeInfo<float>;
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/matrix_free/shape_info.inst.in b/source/matrix_free/shape_info.inst.in
new file mode 100644 (file)
index 0000000..f2c859d
--- /dev/null
@@ -0,0 +1,42 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+for (deal_II_dimension : DIMENSIONS; deal_II_scalar : REAL_SCALARS)
+  {
+    template void internal::MatrixFreeFunctions::ShapeInfo<deal_II_scalar>::
+      reinit<deal_II_dimension>(
+        const Quadrature<1> &,
+        const FiniteElement<deal_II_dimension, deal_II_dimension> &,
+        const unsigned int);
+  }
+
+
+
+for (deal_II_dimension : DIMENSIONS;
+     deal_II_scalar_vectorized : REAL_SCALARS_VECTORIZED)
+  {
+    template void internal::MatrixFreeFunctions::
+      ShapeInfo<deal_II_scalar_vectorized>::reinit<deal_II_dimension>(
+        const Quadrature<1> &,
+        const FiniteElement<deal_II_dimension, deal_II_dimension> &,
+        const unsigned int);
+  }
+
+for (deal_II_scalar_vectorized : REAL_SCALARS_VECTORIZED)
+  {
+    template struct internal::MatrixFreeFunctions::ShapeInfo<
+      deal_II_scalar_vectorized>;
+  }

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.