From: Marc Fehling Date: Mon, 27 Apr 2020 14:09:00 +0000 (+0200) Subject: Precalculate hp::FEValues. X-Git-Tag: v9.2.0-rc1~161^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F9968%2Fhead;p=dealii.git Precalculate hp::FEValues. --- diff --git a/doc/news/changes/minor/20200427Fehling b/doc/news/changes/minor/20200427Fehling new file mode 100644 index 0000000000..dfb04c31c8 --- /dev/null +++ b/doc/news/changes/minor/20200427Fehling @@ -0,0 +1,4 @@ +New: All FEValues objects in a hp::FEValues collection can now be +pre-computed with hp::FEValues::precalculate_fe_values(). +
+(Marc Fehling, 2020/04/27) diff --git a/include/deal.II/hp/fe_values.h b/include/deal.II/hp/fe_values.h index 44b56c85e7..10ffea6706 100644 --- a/include/deal.II/hp/fe_values.h +++ b/include/deal.II/hp/fe_values.h @@ -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 &fe_collection, const QCollection & 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_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 &fe_indices, + const std::vector &mapping_indices, + const std::vector &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> 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 (spacedim != dim ). diff --git a/source/hp/fe_values.cc b/source/hp/fe_values.cc index e980cbe931..ff1c635a00 100644 --- a/source/hp/fe_values.cc +++ b/source/hp/fe_values.cc @@ -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 + #include #include @@ -92,6 +94,64 @@ namespace hp // now there definitely is one! return *fe_values_table(present_fe_values_index); } + + + + template + void + FEValuesBase::precalculate_fe_values( + const std::vector &fe_indices, + const std::vector &mapping_indices, + const std::vector &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( + (*mapping_collection)[mapping_index], + (*fe_collection)[fe_index], + q_collection[q_index], + update_flags); + }); + } + + task_group.join_all(); + } + + + + template + void + FEValuesBase::precalculate_fe_values() + { + const unsigned int size = fe_collection->size(); + std::vector indices(size); + std::iota(indices.begin(), indices.end(), 0); + + precalculate_fe_values(/*fe_indices=*/indices, + /*mapping_indices=*/ + (mapping_collection->size() > 1) ? + indices : + std::vector(size, 0), + /*q_indices=*/ + (q_collection.size() > 1) ? + indices : + std::vector(size, 0)); + } } // namespace hp diff --git a/tests/hp/fe_values_01.cc b/tests/hp/fe_values_01.cc new file mode 100644 index 0000000000..f09f7714f0 --- /dev/null +++ b/tests/hp/fe_values_01.cc @@ -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 + +#include +#include + +#include +#include +#include +#include + +#include "../tests.h" + + + +template +void +test() +{ + const unsigned int max_degree = 3; + + hp::FECollection fes; + for (unsigned int degree = 1; degree <= max_degree; ++degree) + fes.push_back(FE_Q(degree)); + + // provide quadrature and mapping collections corresponding to the size of the + // finite element collection + { + hp::MappingCollection mappings; + hp::QCollection quadratures; + + for (unsigned int i = 0; i < fes.size(); ++i) + { + mappings.push_back(MappingQ(fes[i].degree)); + quadratures.push_back(QGauss(fes[i].degree + 1)); + } + + hp::FEValues 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 quadratures(QGauss(max_degree + 1)); + + hp::FEValues 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 index 0000000000..f445833684 --- /dev/null +++ b/tests/hp/fe_values_01.output @@ -0,0 +1,2 @@ + +DEAL:2d::OK