From 0e2ee3386a0a9f88c2a3b6e48cc6fe46296d2d14 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Wed, 29 Apr 2020 19:13:58 +0200 Subject: [PATCH] Added SmoothnessEstimator::default_fe_series() functions returning default FESeries objects. Removed default constructor and initialize() functions for FESeries objects. --- .../{20200121Fehling-2 => 20200121Fehling} | 0 doc/news/changes/minor/20200121Fehling-1 | 4 -- examples/step-27/step-27.cc | 19 +++---- include/deal.II/fe/fe_series.h | 46 ++------------- .../deal.II/numerics/smoothness_estimator.h | 56 ++++++------------ source/fe/fe_series_fourier.cc | 45 ++------------- source/fe/fe_series_legendre.cc | 43 ++------------ source/numerics/smoothness_estimator.cc | 57 +++++++------------ source/numerics/smoothness_estimator.inst.in | 22 +++---- tests/hp/laplace_mitchell2014_04_peak.cc | 17 ++---- tests/hp/step-27.cc | 13 ++--- tests/mpi/petsc_step-27.cc | 11 ++-- tests/mpi/trilinos_step-27.cc | 11 ++-- tests/numerics/smoothness_estimator_01.cc | 17 ++---- tests/numerics/smoothness_estimator_02.cc | 17 ++---- 15 files changed, 98 insertions(+), 280 deletions(-) rename doc/news/changes/minor/{20200121Fehling-2 => 20200121Fehling} (100%) delete mode 100644 doc/news/changes/minor/20200121Fehling-1 diff --git a/doc/news/changes/minor/20200121Fehling-2 b/doc/news/changes/minor/20200121Fehling similarity index 100% rename from doc/news/changes/minor/20200121Fehling-2 rename to doc/news/changes/minor/20200121Fehling diff --git a/doc/news/changes/minor/20200121Fehling-1 b/doc/news/changes/minor/20200121Fehling-1 deleted file mode 100644 index 15f4f7d57d..0000000000 --- a/doc/news/changes/minor/20200121Fehling-1 +++ /dev/null @@ -1,4 +0,0 @@ -New: Both classes FESeries::Fourier and FESeries::Legendre now offer a -default constructor and can be initialized after their instantiation. -
-(Marc Fehling, 2020/01/21) diff --git a/examples/step-27/step-27.cc b/examples/step-27/step-27.cc index 70d163efdb..87f3bb252d 100644 --- a/examples/step-27/step-27.cc +++ b/examples/step-27/step-27.cc @@ -110,10 +110,10 @@ namespace Step27 hp::QCollection quadrature_collection; hp::QCollection face_quadrature_collection; - hp::QCollection fourier_q_collection; - FESeries::Fourier fourier; - std::vector ln_k; - Table> fourier_coefficients; + hp::QCollection fourier_q_collection; + std::unique_ptr> fourier; + std::vector ln_k; + Table> fourier_coefficients; AffineConstraints constraints; @@ -232,9 +232,8 @@ namespace Step27 // Now we are ready to set-up the FESeries::Fourier object const std::vector n_coefficients_per_direction( fe_collection.size(), N); - fourier.initialize(n_coefficients_per_direction, - fe_collection, - fourier_q_collection); + fourier = std_cxx14::make_unique>( + n_coefficients_per_direction, fe_collection, fourier_q_collection); // We need to resize the matrix of fourier coefficients according to the // number of modes N. @@ -679,9 +678,9 @@ namespace Step27 local_dof_values.reinit(cell->get_fe().dofs_per_cell); cell->get_dof_values(solution, local_dof_values); - fourier.calculate(local_dof_values, - cell->active_fe_index(), - fourier_coefficients); + fourier->calculate(local_dof_values, + cell->active_fe_index(), + fourier_coefficients); // The next thing, as explained in the introduction, is that we wanted // to only fit our exponential decay of Fourier coefficients to the diff --git a/include/deal.II/fe/fe_series.h b/include/deal.II/fe/fe_series.h index dc4fe34234..6a76cf24d8 100644 --- a/include/deal.II/fe/fe_series.h +++ b/include/deal.II/fe/fe_series.h @@ -93,14 +93,6 @@ namespace FESeries public: using CoefficientType = typename std::complex; - /** - * Default constructor. - * - * Leaves all members in an empty or uninitialized state. Please call - * initialize() once you have all data structures ready. - */ - Fourier() = default; - /** * Constructor that initializes all required data structures. * @@ -127,17 +119,6 @@ namespace FESeries const hp::FECollection &fe_collection, const hp::QCollection & q_collection); - /** - * Initialize and overwrite all mandatory data structures for the - * calculation calculation of transformation matrices. - * - * All previously calculated transformation matrices will be cleared. - */ - void - initialize(const std::vector &n_coefficients_per_direction, - const hp::FECollection &fe_collection, - const hp::QCollection & q_collection); - /** * Calculate @p fourier_coefficients of the cell vector field given by * @p local_dof_values corresponding to FiniteElement with @@ -201,7 +182,7 @@ namespace FESeries * Number of coefficients in each direction for each finite element in the * registered hp::FECollection. */ - std::vector n_coefficients_per_direction; + const std::vector n_coefficients_per_direction; /** * hp::FECollection for which transformation matrices will be calculated. @@ -211,7 +192,7 @@ namespace FESeries /** * hp::QCollection used in calculation of transformation matrices. */ - SmartPointer> q_collection; + const hp::QCollection q_collection; /** * Angular frequencies $ 2 \pi {\bf k} $ . @@ -280,14 +261,6 @@ namespace FESeries public: using CoefficientType = double; - /** - * Default constructor. - * - * Leaves all members in an empty or uninitialized state. Please call - * initialize() once you have all data structures ready. - */ - Legendre() = default; - /** * Constructor that initializes all required data structures. * @@ -313,17 +286,6 @@ namespace FESeries const hp::FECollection &fe_collection, const hp::QCollection & q_collection); - /** - * Initialize and overwrite all mandatory data structures for the - * calculation of transformation matrices. - * - * All previously calculated transformation matrices will be cleared. - */ - void - initialize(const std::vector &n_coefficients_per_direction, - const hp::FECollection &fe_collection, - const hp::QCollection & q_collection); - /** * Calculate @p legendre_coefficients of the cell vector field given by * @p local_dof_values corresponding to FiniteElement with @@ -387,7 +349,7 @@ namespace FESeries * Number of coefficients in each direction for each finite element in the * registered hp::FECollection. */ - std::vector n_coefficients_per_direction; + const std::vector n_coefficients_per_direction; /** * hp::FECollection for which transformation matrices will be calculated. @@ -397,7 +359,7 @@ namespace FESeries /** * hp::QCollection used in calculation of transformation matrices. */ - SmartPointer> q_collection; + const hp::QCollection q_collection; /** * Transformation matrices for each FiniteElement. diff --git a/include/deal.II/numerics/smoothness_estimator.h b/include/deal.II/numerics/smoothness_estimator.h index 659ad33049..e73762196d 100644 --- a/include/deal.II/numerics/smoothness_estimator.h +++ b/include/deal.II/numerics/smoothness_estimator.h @@ -225,29 +225,18 @@ namespace SmoothnessEstimator const bool only_flagged_cells = false); /** - * Number of modes for the default configuration of the Legendre series - * expansion FESeries::Legendre. + * Returns a FESeries::Legendre object for Legendre series expansions with + * the default configuration for smoothness estimation purposes. * - * We use as many modes as the highest polynomial of all finite elements - * used plus one, since we start with the first Legendre polynomial which is - * just a constant. - */ - template - std::vector - default_number_of_coefficients_per_direction( - const hp::FECollection &fe_collection); - - /** - * Quadrature collection for the default configuration of the Legendre - * series expansion FESeries::Legendre. - * - * We use Gaussian quadrature designed to yield exact results for the + * For each finite element of the provided @p fe_collection, we use as many + * modes as its polynomial degree plus one, since we start with the first + * Legendre polynomial which is just a constant. Further for each element, + * we use a Gaussian quadrature designed to yield exact results for the * highest order Legendre polynomial used. */ template - hp::QCollection - default_quadrature_collection( - const hp::FECollection &fe_collection); + FESeries::Legendre + default_fe_series(const hp::FECollection &fe_collection); } // namespace Legendre @@ -465,29 +454,18 @@ namespace SmoothnessEstimator const bool only_flagged_cells = false); /** - * Number of modes for the default configuration of the Fourier series - * expansion FESeries::Fourier. - * - * We use as many modes as the highest polynomial degree of all finite - * elements used plus one, and at least three modes. - */ - template - std::vector - default_number_of_coefficients_per_direction( - const hp::FECollection &fe_collection); - - /** - * Quadrature collection for the default configuration of the Fourier series - * expansion FESeries::Fourier. + * Returns a FESeries::Fourier object for Fourier series expansions with + * the default configuration for smoothness estimation purposes. * - * We use a 4-point Gaussian quarature iterated in each dimension by the - * maximal wave number, which is the number of modes decresed by one since - * we start with $k = 0$. + * For each finite element of the provided @p fe_collection, we use as many + * modes as its polynomial degree plus one, and at least three modes. + * Further for each element, we use a 4-point Gaussian quarature iterated in + * each dimension by the maximal wave number, which is the number of modes + * decreased by one since we start with $k = 0$. */ template - hp::QCollection - default_quadrature_collection( - const hp::FECollection &fe_collection); + FESeries::Fourier + default_fe_series(const hp::FECollection &fe_collection); } // namespace Fourier } // namespace SmoothnessEstimator diff --git a/source/fe/fe_series_fourier.cc b/source/fe/fe_series_fourier.cc index e810baa426..814592253b 100644 --- a/source/fe/fe_series_fourier.cc +++ b/source/fe/fe_series_fourier.cc @@ -179,7 +179,7 @@ namespace FESeries const hp::QCollection & q_collection) : n_coefficients_per_direction(n_coefficients_per_direction) , fe_collection(&fe_collection) - , q_collection(&q_collection) + , q_collection(q_collection) , fourier_transform_matrices(fe_collection.size()) { Assert(n_coefficients_per_direction.size() == fe_collection.size() && @@ -219,7 +219,7 @@ namespace FESeries return ( (n_coefficients_per_direction == fourier.n_coefficients_per_direction) && (*fe_collection == *(fourier.fe_collection)) && - (*q_collection == *(fourier.q_collection)) && + (q_collection == fourier.q_collection) && (k_vectors == fourier.k_vectors) && (fourier_transform_matrices == fourier.fourier_transform_matrices)); } @@ -235,7 +235,7 @@ namespace FESeries task_group += Threads::new_task([&, fe]() { ensure_existence(n_coefficients_per_direction, *fe_collection, - *q_collection, + q_collection, k_vectors, fe, fourier_transform_matrices); @@ -246,41 +246,6 @@ namespace FESeries - template - void - Fourier::initialize( - const std::vector & n_coefficients_per_direction, - const hp::FECollection &fe_collection, - const hp::QCollection & q_collection) - { - Assert(n_coefficients_per_direction.size() == fe_collection.size() && - n_coefficients_per_direction.size() == q_collection.size(), - ExcMessage("All parameters are supposed to have the same size.")); - - // set members - this->n_coefficients_per_direction = n_coefficients_per_direction; - this->fe_collection = - SmartPointer>(&fe_collection); - this->q_collection = - SmartPointer>(&q_collection); - - // clean up auxiliary members - k_vectors.reset_values(); - const unsigned int max_n_coefficients_per_direction = - *std::max_element(n_coefficients_per_direction.cbegin(), - n_coefficients_per_direction.cend()); - set_k_vectors(k_vectors, max_n_coefficients_per_direction); - - fourier_transform_matrices.clear(); - fourier_transform_matrices.resize(fe_collection.size()); - - unrolled_coefficients.clear(); - // reserve sufficient memory - unrolled_coefficients.reserve(k_vectors.n_elements()); - } - - - template unsigned int Fourier::get_n_coefficients_per_direction( @@ -299,15 +264,13 @@ namespace FESeries const unsigned int cell_active_fe_index, Table &fourier_coefficients) { - Assert(fe_collection != nullptr && q_collection != nullptr, - ExcMessage("Initialize this FESeries object first!")); for (unsigned int d = 0; d < dim; ++d) AssertDimension(fourier_coefficients.size(d), n_coefficients_per_direction[cell_active_fe_index]); ensure_existence(n_coefficients_per_direction, *fe_collection, - *q_collection, + q_collection, k_vectors, cell_active_fe_index, fourier_transform_matrices); diff --git a/source/fe/fe_series_legendre.cc b/source/fe/fe_series_legendre.cc index 6be7a8a4f0..17e1a7739b 100644 --- a/source/fe/fe_series_legendre.cc +++ b/source/fe/fe_series_legendre.cc @@ -201,7 +201,7 @@ namespace FESeries const hp::QCollection & q_collection) : n_coefficients_per_direction(n_coefficients_per_direction) , fe_collection(&fe_collection) - , q_collection(&q_collection) + , q_collection(q_collection) , legendre_transform_matrices(fe_collection.size()) { Assert(n_coefficients_per_direction.size() == fe_collection.size() && @@ -240,7 +240,7 @@ namespace FESeries return ( (n_coefficients_per_direction == legendre.n_coefficients_per_direction) && (*fe_collection == *(legendre.fe_collection)) && - (*q_collection == *(legendre.q_collection)) && + (q_collection == legendre.q_collection) && (legendre_transform_matrices == legendre.legendre_transform_matrices)); } @@ -255,7 +255,7 @@ namespace FESeries task_group += Threads::new_task([&, fe]() { ensure_existence(n_coefficients_per_direction, *fe_collection, - *q_collection, + q_collection, fe, legendre_transform_matrices); }); @@ -265,39 +265,6 @@ namespace FESeries - template - void - Legendre::initialize( - const std::vector & n_coefficients_per_direction, - const hp::FECollection &fe_collection, - const hp::QCollection & q_collection) - { - Assert(n_coefficients_per_direction.size() == fe_collection.size() && - n_coefficients_per_direction.size() == q_collection.size(), - ExcMessage("All parameters are supposed to have the same size.")); - - // set members - this->n_coefficients_per_direction = n_coefficients_per_direction; - this->fe_collection = - SmartPointer>(&fe_collection); - this->q_collection = - SmartPointer>(&q_collection); - - // clean up auxiliary members - legendre_transform_matrices.clear(); - legendre_transform_matrices.resize(fe_collection.size()); - - unrolled_coefficients.clear(); - // reserve sufficient memory - const unsigned int max_n_coefficients_per_direction = - *std::max_element(n_coefficients_per_direction.cbegin(), - n_coefficients_per_direction.cend()); - unrolled_coefficients.reserve( - Utilities::fixed_power(max_n_coefficients_per_direction)); - } - - - template unsigned int Legendre::get_n_coefficients_per_direction( @@ -316,15 +283,13 @@ namespace FESeries const unsigned int cell_active_fe_index, Table & legendre_coefficients) { - Assert(fe_collection != nullptr && q_collection != nullptr, - ExcMessage("Initialize this FESeries object first!")); for (unsigned int d = 0; d < dim; ++d) AssertDimension(legendre_coefficients.size(d), n_coefficients_per_direction[cell_active_fe_index]); ensure_existence(n_coefficients_per_direction, *fe_collection, - *q_collection, + q_collection, cell_active_fe_index, legendre_transform_matrices); diff --git a/source/numerics/smoothness_estimator.cc b/source/numerics/smoothness_estimator.cc index 568765d245..9a22348534 100644 --- a/source/numerics/smoothness_estimator.cc +++ b/source/numerics/smoothness_estimator.cc @@ -286,27 +286,16 @@ namespace SmoothnessEstimator template - std::vector - default_number_of_coefficients_per_direction( - const hp::FECollection &fe_collection) + FESeries::Legendre + default_fe_series(const hp::FECollection &fe_collection) { + // Default number of coefficients per direction. std::vector n_coefficients_per_direction; for (unsigned int i = 0; i < fe_collection.size(); ++i) n_coefficients_per_direction.push_back(fe_collection[i].degree + 1); - return n_coefficients_per_direction; - } - - - - template - hp::QCollection - default_quadrature_collection( - const hp::FECollection &fe_collection) - { - const std::vector n_modes = - default_number_of_coefficients_per_direction(fe_collection); - + // Default quadrature collection. + // // We initialize a FESeries::Legendre expansion object object which will // be used to calculate the expansion coefficients. In addition to the // hp::FECollection, we need to provide quadrature rules hp::QCollection @@ -317,18 +306,20 @@ namespace SmoothnessEstimator // rule. As a default, we use the same quadrature formula for each finite // element, namely a Gauss formula that yields exact results for the // highest order Legendre polynomial used. - + // // We start with the zeroth Legendre polynomial which is just a constant, // so the highest Legendre polynomial will be of order (n_modes - 1). hp::QCollection q_collection; for (unsigned int i = 0; i < fe_collection.size(); ++i) { - const QGauss quadrature(n_modes[i]); + const QGauss quadrature(n_coefficients_per_direction[i]); const QSorted quadrature_sorted(quadrature); q_collection.push_back(quadrature_sorted); } - return q_collection; + return FESeries::Legendre(n_coefficients_per_direction, + fe_collection, + q_collection); } } // namespace Legendre @@ -578,28 +569,17 @@ namespace SmoothnessEstimator template - std::vector - default_number_of_coefficients_per_direction( - const hp::FECollection &fe_collection) + FESeries::Fourier + default_fe_series(const hp::FECollection &fe_collection) { + // Default number of coefficients per direction. std::vector n_coefficients_per_direction; for (unsigned int i = 0; i < fe_collection.size(); ++i) n_coefficients_per_direction.push_back( std::max(3, fe_collection[i].degree + 1)); - return n_coefficients_per_direction; - } - - - - template - hp::QCollection - default_quadrature_collection( - const hp::FECollection &fe_collection) - { - const std::vector n_modes = - default_number_of_coefficients_per_direction(fe_collection); - + // Default quadrature collection. + // // We initialize a series expansion object object which will be used to // calculate the expansion coefficients. In addition to the // hp::FECollection, we need to provide quadrature rules hp::QCollection @@ -616,12 +596,15 @@ namespace SmoothnessEstimator hp::QCollection q_collection; for (unsigned int i = 0; i < fe_collection.size(); ++i) { - const QIterated quadrature(base_quadrature, n_modes[i] - 1); + const QIterated quadrature(base_quadrature, + n_coefficients_per_direction[i] - 1); const QSorted quadrature_sorted(quadrature); q_collection.push_back(quadrature_sorted); } - return q_collection; + return FESeries::Fourier(n_coefficients_per_direction, + fe_collection, + q_collection); } } // namespace Fourier } // namespace SmoothnessEstimator diff --git a/source/numerics/smoothness_estimator.inst.in b/source/numerics/smoothness_estimator.inst.in index 912fa2d70c..030466be4c 100644 --- a/source/numerics/smoothness_estimator.inst.in +++ b/source/numerics/smoothness_estimator.inst.in @@ -64,20 +64,14 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { #if deal_II_dimension <= deal_II_space_dimension - template std::vector SmoothnessEstimator::Legendre:: - default_number_of_coefficients_per_direction( - const hp::FECollection &); - template hp::QCollection SmoothnessEstimator::Legendre:: - default_quadrature_collection( - const hp::FECollection &); + template FESeries::Legendre + SmoothnessEstimator::Legendre::default_fe_series( + const hp::FECollection &); - template std::vector SmoothnessEstimator::Fourier:: - default_number_of_coefficients_per_direction( - const hp::FECollection &); - template hp::QCollection SmoothnessEstimator::Fourier:: - default_quadrature_collection( - const hp::FECollection &); + template FESeries::Fourier + SmoothnessEstimator::Fourier::default_fe_series( + const hp::FECollection &); #endif } diff --git a/tests/hp/laplace_mitchell2014_04_peak.cc b/tests/hp/laplace_mitchell2014_04_peak.cc index 26d7c49527..2b16cfdd8c 100644 --- a/tests/hp/laplace_mitchell2014_04_peak.cc +++ b/tests/hp/laplace_mitchell2014_04_peak.cc @@ -111,9 +111,7 @@ private: hp::QCollection quadrature_face; - std::vector n_coefficients_per_direction; - hp::QCollection expansion_q_collection; - FESeries::Legendre legendre; + std::unique_ptr> legendre; }; template @@ -143,15 +141,8 @@ Problem4::Problem4(const Function &force_function, // after the FECollection has been generated, create a corresponding legendre // series expansion object - n_coefficients_per_direction = - SmoothnessEstimator::Legendre::default_number_of_coefficients_per_direction( - Laplace::fe); - expansion_q_collection = - SmoothnessEstimator::Legendre::default_quadrature_collection( - Laplace::fe); - legendre.initialize(n_coefficients_per_direction, - Laplace::fe, - expansion_q_collection); + legendre = std_cxx14::make_unique>( + SmoothnessEstimator::Legendre::default_fe_series(Laplace::fe)); } @@ -163,7 +154,7 @@ Problem4::substitute_h_for_p() Vector smoothness_indicators( Laplace::triangulation.n_active_cells()); SmoothnessEstimator::Legendre::coefficient_decay_per_direction( - legendre, + *legendre, Laplace::dof_handler, Laplace::solution, smoothness_indicators); diff --git a/tests/hp/step-27.cc b/tests/hp/step-27.cc index eafadd5c0e..6fbe477f99 100644 --- a/tests/hp/step-27.cc +++ b/tests/hp/step-27.cc @@ -104,9 +104,9 @@ namespace Step27 Vector solution; Vector system_rhs; - const unsigned int max_degree; - hp::QCollection fourier_q_collection; - FESeries::Fourier fourier; + const unsigned int max_degree; + hp::QCollection fourier_q_collection; + std::unique_ptr> fourier; }; @@ -155,9 +155,8 @@ namespace Step27 const std::vector n_coefficients_per_direction( fe_collection.size(), max_degree); - fourier.initialize(n_coefficients_per_direction, - fe_collection, - fourier_q_collection); + fourier = std_cxx14::make_unique>( + n_coefficients_per_direction, fe_collection, fourier_q_collection); } @@ -292,7 +291,7 @@ namespace Step27 Vector smoothness_indicators; SmoothnessEstimator::Fourier::coefficient_decay( - fourier, + *fourier, dof_handler, solution, smoothness_indicators, diff --git a/tests/mpi/petsc_step-27.cc b/tests/mpi/petsc_step-27.cc index eb95c3caff..f461b06a80 100644 --- a/tests/mpi/petsc_step-27.cc +++ b/tests/mpi/petsc_step-27.cc @@ -107,8 +107,8 @@ namespace Step27 hp::QCollection quadrature_collection; hp::QCollection face_quadrature_collection; - hp::QCollection fourier_q_collection; - FESeries::Fourier fourier; + hp::QCollection fourier_q_collection; + std::unique_ptr> fourier; AffineConstraints constraints; @@ -188,9 +188,8 @@ namespace Step27 const std::vector n_coefficients_per_direction( fe_collection.size(), max_degree); - fourier.initialize(n_coefficients_per_direction, - fe_collection, - fourier_q_collection); + fourier = std_cxx14::make_unique>( + n_coefficients_per_direction, fe_collection, fourier_q_collection); } @@ -376,7 +375,7 @@ namespace Step27 Vector smoothness_indicators(triangulation.n_active_cells()); SmoothnessEstimator::Fourier::coefficient_decay( - fourier, + *fourier, dof_handler, solution, smoothness_indicators, diff --git a/tests/mpi/trilinos_step-27.cc b/tests/mpi/trilinos_step-27.cc index 692a3c2c61..555936712a 100644 --- a/tests/mpi/trilinos_step-27.cc +++ b/tests/mpi/trilinos_step-27.cc @@ -107,8 +107,8 @@ namespace Step27 hp::QCollection quadrature_collection; hp::QCollection face_quadrature_collection; - hp::QCollection fourier_q_collection; - FESeries::Fourier fourier; + hp::QCollection fourier_q_collection; + std::unique_ptr> fourier; AffineConstraints constraints; @@ -188,9 +188,8 @@ namespace Step27 const std::vector n_coefficients_per_direction( fe_collection.size(), max_degree); - fourier.initialize(n_coefficients_per_direction, - fe_collection, - fourier_q_collection); + fourier = std_cxx14::make_unique>( + n_coefficients_per_direction, fe_collection, fourier_q_collection); } @@ -377,7 +376,7 @@ namespace Step27 Vector smoothness_indicators(triangulation.n_active_cells()); SmoothnessEstimator::Fourier::coefficient_decay( - fourier, + *fourier, dof_handler, solution, smoothness_indicators, diff --git a/tests/numerics/smoothness_estimator_01.cc b/tests/numerics/smoothness_estimator_01.cc index a08429ace5..f5bd97c58d 100644 --- a/tests/numerics/smoothness_estimator_01.cc +++ b/tests/numerics/smoothness_estimator_01.cc @@ -179,11 +179,12 @@ test(const LegendreFunction &func, const unsigned int poly_degree) for (unsigned int p = 1; p <= max_poly; p++) fe_collection.push_back(FE_Q(p)); - const unsigned int fe_index = poly_degree - 1; - const std::vector n_coefficients_per_direction = - SmoothnessEstimator::Legendre::default_number_of_coefficients_per_direction( - fe_collection); - const unsigned int n_modes = n_coefficients_per_direction[fe_index]; + FESeries::Legendre legendre = + SmoothnessEstimator::Legendre::default_fe_series(fe_collection); + + const unsigned int fe_index = poly_degree - 1; + const unsigned int n_modes = + legendre.get_n_coefficients_per_direction(fe_index); // custom predicate: // p-ref for linear elements and use j=1,...,pe otherwise. @@ -201,12 +202,6 @@ test(const LegendreFunction &func, const unsigned int poly_degree) Vector values(dof_handler.n_dofs()); VectorTools::interpolate(dof_handler, func, values); - hp::QCollection q_collection = - SmoothnessEstimator::Legendre::default_quadrature_collection(fe_collection); - FESeries::Legendre legendre(n_coefficients_per_direction, - fe_collection, - q_collection); - const Table &coeff_in = func.get_coefficients(); Table coeff_out; resize(coeff_out, n_modes); diff --git a/tests/numerics/smoothness_estimator_02.cc b/tests/numerics/smoothness_estimator_02.cc index 0b59a5299b..6c007b521a 100644 --- a/tests/numerics/smoothness_estimator_02.cc +++ b/tests/numerics/smoothness_estimator_02.cc @@ -193,11 +193,12 @@ test(const unsigned int poly_degree) for (unsigned int p = 1; p <= max_poly; ++p) fe_collection.push_back(FE_Q(p)); - const unsigned int fe_index = poly_degree - 1; - const std::vector n_coefficients_per_direction = - SmoothnessEstimator::Fourier::default_number_of_coefficients_per_direction( - fe_collection); - const unsigned int n_modes = n_coefficients_per_direction[fe_index]; + FESeries::Fourier fourier = + SmoothnessEstimator::Fourier::default_fe_series(fe_collection); + + const unsigned int fe_index = poly_degree - 1; + const unsigned int n_modes = + fourier.get_n_coefficients_per_direction(fe_index); Assert((poly_degree >= 1) && (poly_degree <= max_poly), ExcInternalError()); Assert((n_modes >= 3) && (n_modes <= max_poly + 1), ExcInternalError()); @@ -277,12 +278,6 @@ test(const unsigned int poly_degree) local_dof_values.reinit(cell->get_fe().dofs_per_cell); cell->get_dof_values(values, local_dof_values); - hp::QCollection q_collection = - SmoothnessEstimator::Fourier::default_quadrature_collection(fe_collection); - FESeries::Fourier fourier(n_coefficients_per_direction, - fe_collection, - q_collection); - Table> coeff_out; coeff_out.reinit(size); fourier.calculate(local_dof_values, cell->active_fe_index(), coeff_out); -- 2.39.5