]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Precalculate hp::FEValues. 9968/head
authorMarc Fehling <marc.fehling@gmx.net>
Mon, 27 Apr 2020 14:09:00 +0000 (16:09 +0200)
committerMarc Fehling <marc.fehling@gmx.net>
Wed, 29 Apr 2020 08:32:10 +0000 (10:32 +0200)
doc/news/changes/minor/20200427Fehling [new file with mode: 0644]
include/deal.II/hp/fe_values.h
source/hp/fe_values.cc
tests/hp/fe_values_01.cc [new file with mode: 0644]
tests/hp/fe_values_01.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20200427Fehling b/doc/news/changes/minor/20200427Fehling
new file mode 100644 (file)
index 0000000..dfb04c3
--- /dev/null
@@ -0,0 +1,4 @@
+New: All FEValues objects in a hp::FEValues collection can now be
+pre-computed with hp::FEValues::precalculate_fe_values().
+<br>
+(Marc Fehling, 2020/04/27)
index 44b56c85e79484af020c53b22803befad1fa603d..10ffea670633de12ab973537e02e5a4352b868ec 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2003 - 2018 by the deal.II authors
+// Copyright (C) 2003 - 2020 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -44,9 +44,12 @@ namespace hp
    * that is common to them. The main task of this class is to provide a
    * table where for every combination of finite element, mapping, and
    * quadrature object from their corresponding collection objects there is
-   * a matching ::FEValues, ::FEFaceValues, or ::FESubfaceValues object. To
-   * make things more efficient, however, these FE*Values objects are only
-   * created once requested (lazy allocation).
+   * a matching ::FEValues, ::FEFaceValues, or ::FESubfaceValues object.
+   *
+   * To make things more efficient, however, these FE*Values objects are only
+   * created once requested (lazy allocation). Alternatively if desired, this
+   * can be bypassed by computing all objects in advance with the corresponding
+   * precalculate_fe_values() function.
    *
    * The first template parameter denotes the space dimension we are in, the
    * second the dimensionality of the object that we integrate on, i.e. for
@@ -73,6 +76,7 @@ namespace hp
       const FECollection<dim, FEValuesType::space_dimension> &fe_collection,
       const QCollection<q_dim> &                              q_collection,
       const UpdateFlags                                       update_flags);
+
     /**
      * Constructor. This constructor is equivalent to the other one except
      * that it makes the object use a $Q_1$ mapping (i.e., an object of type
@@ -83,6 +87,35 @@ namespace hp
       const QCollection<q_dim> &                              q_collection,
       const UpdateFlags                                       update_flags);
 
+    /**
+     * For timing purposes it may be useful to create all required FE*Values
+     * objects in advance, rather than computing them on request via lazy
+     * allocation as usual in this class.
+     *
+     * This function precalculates the FE*Values objects corresponding to the
+     * provided parameters: The total of all vector entries corresponding to the
+     * same index describes an FE*Values object similarly to select_fe_values().
+     */
+    void
+    precalculate_fe_values(const std::vector<unsigned int> &fe_indices,
+                           const std::vector<unsigned int> &mapping_indices,
+                           const std::vector<unsigned int> &q_indices);
+
+    /**
+     * Same as above, geared to the most common use of hp::FEValues objects in
+     * which FE, quadrature and mapping indices are similar on each individual
+     * cell.
+     *
+     * FE*Values objects are created for every FE in the FECollection, with
+     * quadrature and mapping corresponding to the same index from the
+     * QuadratureCollection and MappingCollection, respectively.
+     *
+     * If QuadratureCollection or MappingCollection contains only one object, it
+     * is used for all FE*Values objects.
+     */
+    void
+    precalculate_fe_values();
+
     /**
      * Get a reference to the collection of finite element objects used
      * here.
@@ -161,7 +194,7 @@ namespace hp
      * within the q_collection.
      *
      * Initially, all entries have zero pointers, and we will allocate them
-     * lazily as needed in select_fe_values().
+     * lazily as needed in select_fe_values() or precalculate_fe_values().
      */
     Table<3, std::shared_ptr<FEValuesType>> fe_values_table;
 
@@ -222,7 +255,9 @@ namespace hp
    * Note that ::FEValues objects are created on the fly, i.e. only as they
    * are needed. This ensures that we do not create objects for every
    * combination of finite element, quadrature formula and mapping, but only
-   * those that will actually be needed.
+   * those that will actually be needed. Alternatively if desired, this
+   * can be bypassed by computing all objects in advance with the corresponding
+   * hp::FEValuesBase::precalculate_fe_values() function.
    *
    * This class has not yet been implemented for the use in the codimension
    * one case (<tt>spacedim != dim </tt>).
index e980cbe931bb8690a7113c2032554c5a1a609835..ff1c635a00d75bc85987931741ebc6b460d059f1 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2003 - 2018 by the deal.II authors
+// Copyright (C) 2003 - 2020 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -13,6 +13,8 @@
 //
 // ---------------------------------------------------------------------
 
+#include <deal.II/base/thread_management.h>
+
 #include <deal.II/fe/mapping_q1.h>
 
 #include <deal.II/hp/fe_values.h>
@@ -92,6 +94,64 @@ namespace hp
     // now there definitely is one!
     return *fe_values_table(present_fe_values_index);
   }
+
+
+
+  template <int dim, int q_dim, class FEValuesType>
+  void
+  FEValuesBase<dim, q_dim, FEValuesType>::precalculate_fe_values(
+    const std::vector<unsigned int> &fe_indices,
+    const std::vector<unsigned int> &mapping_indices,
+    const std::vector<unsigned int> &q_indices)
+  {
+    AssertDimension(fe_indices.size(), mapping_indices.size());
+    AssertDimension(fe_indices.size(), q_indices.size());
+
+    Threads::TaskGroup<> task_group;
+    for (unsigned int i = 0; i < fe_indices.size(); ++i)
+      {
+        const unsigned int fe_index      = fe_indices[i],
+                           mapping_index = mapping_indices[i],
+                           q_index       = q_indices[i];
+
+        AssertIndexRange(fe_index, fe_collection->size());
+        AssertIndexRange(mapping_index, mapping_collection->size());
+        AssertIndexRange(q_index, q_collection.size());
+
+        task_group +=
+          Threads::new_task([&, fe_index, mapping_index, q_index]() {
+            fe_values_table(TableIndices<3>(fe_index, mapping_index, q_index)) =
+              std::make_shared<FEValuesType>(
+                (*mapping_collection)[mapping_index],
+                (*fe_collection)[fe_index],
+                q_collection[q_index],
+                update_flags);
+          });
+      }
+
+    task_group.join_all();
+  }
+
+
+
+  template <int dim, int q_dim, class FEValuesType>
+  void
+  FEValuesBase<dim, q_dim, FEValuesType>::precalculate_fe_values()
+  {
+    const unsigned int        size = fe_collection->size();
+    std::vector<unsigned int> indices(size);
+    std::iota(indices.begin(), indices.end(), 0);
+
+    precalculate_fe_values(/*fe_indices=*/indices,
+                           /*mapping_indices=*/
+                           (mapping_collection->size() > 1) ?
+                             indices :
+                             std::vector<unsigned int>(size, 0),
+                           /*q_indices=*/
+                           (q_collection.size() > 1) ?
+                             indices :
+                             std::vector<unsigned int>(size, 0));
+  }
 } // namespace hp
 
 
diff --git a/tests/hp/fe_values_01.cc b/tests/hp/fe_values_01.cc
new file mode 100644 (file)
index 0000000..f09f771
--- /dev/null
@@ -0,0 +1,82 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// precalculate hp::FEValues objects
+
+
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/fe_values.h>
+#include <deal.II/hp/mapping_collection.h>
+#include <deal.II/hp/q_collection.h>
+
+#include "../tests.h"
+
+
+
+template <int dim>
+void
+test()
+{
+  const unsigned int max_degree = 3;
+
+  hp::FECollection<dim> fes;
+  for (unsigned int degree = 1; degree <= max_degree; ++degree)
+    fes.push_back(FE_Q<dim>(degree));
+
+  // provide quadrature and mapping collections corresponding to the size of the
+  // finite element collection
+  {
+    hp::MappingCollection<dim> mappings;
+    hp::QCollection<dim>       quadratures;
+
+    for (unsigned int i = 0; i < fes.size(); ++i)
+      {
+        mappings.push_back(MappingQ<dim>(fes[i].degree));
+        quadratures.push_back(QGauss<dim>(fes[i].degree + 1));
+      }
+
+    hp::FEValues<dim> hp_fe_values(mappings, fes, quadratures, update_values);
+    hp_fe_values.precalculate_fe_values();
+  }
+
+  // provide just one quadrature rule and the default mapping object
+  {
+    hp::QCollection<dim> quadratures(QGauss<dim>(max_degree + 1));
+
+    hp::FEValues<dim> hp_fe_values(fes, quadratures, update_values);
+    hp_fe_values.precalculate_fe_values();
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+int
+main()
+{
+  initlog();
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+}
diff --git a/tests/hp/fe_values_01.output b/tests/hp/fe_values_01.output
new file mode 100644 (file)
index 0000000..f445833
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL:2d::OK

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.