From eacb525b04ede231b493229d19054d66c63bfa10 Mon Sep 17 00:00:00 2001 From: kronbichler Date: Fri, 4 Apr 2014 12:37:11 +0000 Subject: [PATCH] Fix DoFTools::extract_constant_modes for FE_Q_DG0 git-svn-id: https://svn.dealii.org/trunk@32714 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 8 ++ deal.II/include/deal.II/fe/fe.h | 23 ++-- deal.II/include/deal.II/fe/fe_dgp.h | 3 +- deal.II/include/deal.II/fe/fe_dgq.h | 3 +- deal.II/include/deal.II/fe/fe_face.h | 11 +- deal.II/include/deal.II/fe/fe_nedelec.h | 3 +- deal.II/include/deal.II/fe/fe_q_base.h | 3 +- deal.II/include/deal.II/fe/fe_q_dg0.h | 9 ++ .../include/deal.II/fe/fe_q_hierarchical.h | 3 +- .../include/deal.II/fe/fe_raviart_thomas.h | 7 +- deal.II/include/deal.II/fe/fe_system.h | 3 +- deal.II/source/dofs/dof_tools.cc | 47 +++++--- tests/fe/element_constant_modes.cc | 6 +- tests/fe/element_constant_modes.output | 36 ++++++ tests/mpi/extract_constant_modes_02.cc | 104 ++++++++++++++++++ .../extract_constant_modes_02.mpirun=4.output | 4 + 16 files changed, 237 insertions(+), 36 deletions(-) create mode 100644 tests/mpi/extract_constant_modes_02.cc create mode 100644 tests/mpi/extract_constant_modes_02.mpirun=4.output diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index c45afc16b7..a49dcacfe2 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -148,6 +148,14 @@ inconvenience this causes.

Specific improvements

    +
  1. Fixed: DoFTools::extract_constant_modes now correctly identifies both + constant modes in the scalar element FE_Q_DG0, which has been realized by a + few modifications in how the constant modes propagate from the element to + the extract_constant_modes function. +
    + (Martin Kronbichler, 2014/04/04) +
  2. +
  3. Fixed: GridTools::laplace_transform had previously announced in the documentation that one can also set the location of interior points, but this was not in fact what was implemented. This has now been fixed: diff --git a/deal.II/include/deal.II/fe/fe.h b/deal.II/include/deal.II/fe/fe.h index 84604ffc7e..1be926e630 100644 --- a/deal.II/include/deal.II/fe/fe.h +++ b/deal.II/include/deal.II/fe/fe.h @@ -1314,13 +1314,22 @@ public: block_mask (const ComponentMask &component_mask) const; /** - * Returns a list of constant modes of the element. The returns table has as - * many rows as there are components in the element and dofs_per_cell - * columns. To each component of the finite element, the row in the returned - * table contains a basis representation of the constant function 1 on the - * element. - */ - virtual Table<2,bool> get_constant_modes () const; + * Returns a list of constant modes of the element. The number of rows in + * the resulting table depends on the elements in use. For standard + * elements, the table has as many rows as there are components in the + * element and dofs_per_cell columns. To each component of the finite + * element, the row in the returned table contains a basis representation of + * the constant function 1 on the element. However, there are some scalar + * elements where there is more than one constant mode, e.g. the element + * FE_Q_DG0. + * + * In order to match the constant modes to the actual components in the + * element, the returned data structure also returns a vector with as many + * components as there are constant modes on the element that contains the + * component number. + */ + virtual std::pair,std::vector > + get_constant_modes () const; //@} diff --git a/deal.II/include/deal.II/fe/fe_dgp.h b/deal.II/include/deal.II/fe/fe_dgp.h index 6343a7c637..9fd6a5f04a 100644 --- a/deal.II/include/deal.II/fe/fe_dgp.h +++ b/deal.II/include/deal.II/fe/fe_dgp.h @@ -331,7 +331,8 @@ public: * Returns a list of constant modes of the element. For this element, the * first entry is true, all other are false. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; protected: diff --git a/deal.II/include/deal.II/fe/fe_dgq.h b/deal.II/include/deal.II/fe/fe_dgq.h index 103c135c7e..ec0c7f9e18 100644 --- a/deal.II/include/deal.II/fe/fe_dgq.h +++ b/deal.II/include/deal.II/fe/fe_dgq.h @@ -331,7 +331,8 @@ public: * Returns a list of constant modes of the element. For this element, it * simply returns one row with all entries set to true. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; /** * Determine an estimate for the diff --git a/deal.II/include/deal.II/fe/fe_face.h b/deal.II/include/deal.II/fe/fe_face.h index f6bd4606b5..2ee366f50f 100644 --- a/deal.II/include/deal.II/fe/fe_face.h +++ b/deal.II/include/deal.II/fe/fe_face.h @@ -128,7 +128,8 @@ public: * Returns a list of constant modes of the element. For this element, it * simply returns one row with all entries set to true. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; private: /** @@ -230,7 +231,8 @@ public: * Returns a list of constant modes of the element. For this element, it * simply returns one row with all entries set to true. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; protected: virtual @@ -437,9 +439,10 @@ public: * Returns a list of constant modes of the element. For this element, the * first entry on each face is true, all other are false (as the constant * function is represented by the first base function of Legendre - * polynomials. + * polynomials). */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; private: /** diff --git a/deal.II/include/deal.II/fe/fe_nedelec.h b/deal.II/include/deal.II/fe/fe_nedelec.h index 81dad103af..f09e5010bd 100644 --- a/deal.II/include/deal.II/fe/fe_nedelec.h +++ b/deal.II/include/deal.II/fe/fe_nedelec.h @@ -285,7 +285,8 @@ public: /** * Returns a list of constant modes of the element. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; virtual std::size_t memory_consumption () const; virtual FiniteElement *clone() const; diff --git a/deal.II/include/deal.II/fe/fe_q_base.h b/deal.II/include/deal.II/fe/fe_q_base.h index 6870f28fd2..6bfb8a2ace 100644 --- a/deal.II/include/deal.II/fe/fe_q_base.h +++ b/deal.II/include/deal.II/fe/fe_q_base.h @@ -199,7 +199,8 @@ public: * Returns a list of constant modes of the element. For this element, the * list consists of true arguments for all components. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; /** * @name Functions to support hp diff --git a/deal.II/include/deal.II/fe/fe_q_dg0.h b/deal.II/include/deal.II/fe/fe_q_dg0.h index 6f0f339e01..90b28f579a 100644 --- a/deal.II/include/deal.II/fe/fe_q_dg0.h +++ b/deal.II/include/deal.II/fe/fe_q_dg0.h @@ -312,6 +312,15 @@ public: virtual bool has_support_on_face (const unsigned int shape_index, const unsigned int face_index) const; + /** + * Returns a list of constant modes of the element. For this element, there + * are two constant modes despite the element is scalar: The first constant + * mode is all ones for the usual FE_Q basis and the second one only using + * the discontinuous part. + */ + virtual std::pair, std::vector > + get_constant_modes () const; + protected: /** * @p clone function instead of a copy constructor. diff --git a/deal.II/include/deal.II/fe/fe_q_hierarchical.h b/deal.II/include/deal.II/fe/fe_q_hierarchical.h index fd6602c5b1..64d3b644cf 100644 --- a/deal.II/include/deal.II/fe/fe_q_hierarchical.h +++ b/deal.II/include/deal.II/fe/fe_q_hierarchical.h @@ -386,7 +386,8 @@ public: * list consists of true arguments for the first vertex shape functions and * false for the remaining ones. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; protected: /** diff --git a/deal.II/include/deal.II/fe/fe_raviart_thomas.h b/deal.II/include/deal.II/fe/fe_raviart_thomas.h index bda93239da..01182aef01 100644 --- a/deal.II/include/deal.II/fe/fe_raviart_thomas.h +++ b/deal.II/include/deal.II/fe/fe_raviart_thomas.h @@ -148,10 +148,11 @@ public: const VectorSlice > > &values) const; /** - * Returns a list of constant modes of the element. For this element, the - * list consists of true arguments for all components. + * Returns a list of constant modes of the element. This method is currently + * not correctly implemented because it returns ones for all components. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; virtual std::size_t memory_consumption () const; virtual FiniteElement *clone() const; diff --git a/deal.II/include/deal.II/fe/fe_system.h b/deal.II/include/deal.II/fe/fe_system.h index 7022cdd694..04f9736688 100644 --- a/deal.II/include/deal.II/fe/fe_system.h +++ b/deal.II/include/deal.II/fe/fe_system.h @@ -476,7 +476,8 @@ public: * table contains a basis representation of the constant function 1 on the * element. Concatenates the constant modes of each base element. */ - virtual Table<2,bool> get_constant_modes () const; + virtual std::pair, std::vector > + get_constant_modes () const; /** * @name Functions to support hp diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc index 3ea5499334..1ca207eb0d 100644 --- a/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/source/dofs/dof_tools.cc @@ -984,12 +984,6 @@ namespace DoFTools Assert (component_mask.represents_n_components(n_components), ExcDimensionMismatch(n_components, component_mask.size())); - std::vector localized_component (n_components, - numbers::invalid_unsigned_int); - unsigned int n_components_selected = 0; - for (unsigned int i=0; i dofs_by_component (dof_handler.n_locally_owned_dofs()); internal::get_component_association (dof_handler, component_mask, @@ -1008,17 +1002,37 @@ namespace DoFTools if (component_mask[dofs_by_component[i]]) component_numbering[i] = count++; - // First count the number of dofs in the current component. - constant_modes.clear (); - constant_modes.resize (n_components_selected, std::vector(n_selected_dofs, - false)); - - // Loop over all owned cells and ask the element for the constant modes + // get the element constant modes and find a translation table between + // index in the constant modes and the components. + // + // TODO: We might be able to extend this also for elements which do not + // have the same constant modes, but that is messy... const dealii::hp::FECollection fe_collection (dof_handler.get_fe()); std::vector > element_constant_modes; + std::vector > > + constant_mode_to_component_translation(n_components); + unsigned int n_constant_modes = 0; for (unsigned int f=0; f, std::vector > data + = fe_collection[f].get_constant_modes(); + element_constant_modes.push_back(data.first); + if (f==0) + for (unsigned int i=0; i(n_selected_dofs, + false)); + + // Loop over all owned cells and ask the element for the constant modes typename DH::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); @@ -1036,8 +1050,11 @@ namespace DoFTools locally_owned_dofs.index_within_set(dof_indices[i]); const unsigned int comp = dofs_by_component[loc_index]; if (component_mask[comp]) - constant_modes[localized_component[comp]][component_numbering[loc_index]] = - element_constant_modes[cell->active_fe_index()](comp,i); + for (unsigned int j=0; jactive_fe_index()] + (constant_mode_to_component_translation[comp][j].second,i); } } } diff --git a/tests/fe/element_constant_modes.cc b/tests/fe/element_constant_modes.cc index 4ae2eab1fa..54e7d407c3 100644 --- a/tests/fe/element_constant_modes.cc +++ b/tests/fe/element_constant_modes.cc @@ -21,6 +21,7 @@ #include #include #include +#include #include #include #include @@ -30,9 +31,9 @@ template void print_constant_modes(const FiniteElement &fe) { - Table<2,bool> constant_modes = fe.get_constant_modes(); deallog << "Testing " << fe.get_name() << std::endl; + Table<2,bool> constant_modes = fe.get_constant_modes().first; for (unsigned int r=0; r(1)); print_constant_modes(FESystem(FE_Q(1), 2, FE_Q(2), 1)); print_constant_modes(FESystem(FE_DGP(1), 1, FE_Q_iso_Q1(2), 1)); + print_constant_modes(FE_Q_DG0(1)); + print_constant_modes(FESystem(FE_Q_DG0(2), 1, FE_Q(1), 2)); + print_constant_modes(FESystem(FE_Q(1), 2, FE_Q_DG0(1), 2)); } template <> diff --git a/tests/fe/element_constant_modes.output b/tests/fe/element_constant_modes.output index 428dbf5c6e..076019878b 100644 --- a/tests/fe/element_constant_modes.output +++ b/tests/fe/element_constant_modes.output @@ -44,6 +44,24 @@ DEAL::Testing FESystem<2>[FE_DGP<2>(1)-FE_Q_iso_Q1<2>(2)] DEAL::0 0 0 0 0 0 0 0 1 0 0 0 DEAL::1 1 1 1 1 1 1 1 0 0 0 1 DEAL:: +DEAL::Testing FE_Q_DG0<2>(1) +DEAL::1 1 1 1 0 +DEAL::0 0 0 0 1 +DEAL:: +DEAL::Testing FESystem<2>[FE_Q_DG0<2>(2)-FE_Q<2>(1)^2] +DEAL::1 0 0 1 0 0 1 0 0 1 0 0 1 1 1 1 1 0 +DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 +DEAL::0 1 0 0 1 0 0 1 0 0 1 0 0 0 0 0 0 0 +DEAL::0 0 1 0 0 1 0 0 1 0 0 1 0 0 0 0 0 0 +DEAL:: +DEAL::Testing FESystem<2>[FE_Q<2>(1)^2-FE_Q_DG0<2>(1)^2] +DEAL::1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 0 0 +DEAL::0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 0 +DEAL::0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 +DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 +DEAL::0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 +DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 +DEAL:: DEAL::Testing FE_Q<3>(1) DEAL::1 1 1 1 1 1 1 1 DEAL:: @@ -77,3 +95,21 @@ DEAL::Testing FESystem<3>[FE_DGP<3>(1)-FE_Q_iso_Q1<3>(2)] DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 DEAL::1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 1 DEAL:: +DEAL::Testing FE_Q_DG0<3>(1) +DEAL::1 1 1 1 1 1 1 1 0 +DEAL::0 0 0 0 0 0 0 0 1 +DEAL:: +DEAL::Testing FESystem<3>[FE_Q_DG0<3>(2)-FE_Q<3>(1)^2] +DEAL::1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 +DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 +DEAL::0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +DEAL::0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +DEAL:: +DEAL::Testing FESystem<3>[FE_Q<3>(1)^2-FE_Q_DG0<3>(1)^2] +DEAL::1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 0 0 +DEAL::0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 0 +DEAL::0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 +DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 +DEAL::0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0 +DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 +DEAL:: diff --git a/tests/mpi/extract_constant_modes_02.cc b/tests/mpi/extract_constant_modes_02.cc new file mode 100644 index 0000000000..d03cd6002b --- /dev/null +++ b/tests/mpi/extract_constant_modes_02.cc @@ -0,0 +1,104 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2009 - 2013 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + + +// test DoFTools::extract_constant_modes with FE_Q_DG0 + +#include "../tests.h" +#include "coarse_grid_common.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + + +template +void test() +{ + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + std::vector sub(2); + sub[0] = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + sub[1] = 1; + GridGenerator::subdivided_hyper_rectangle(static_cast&>(tr), + sub, Point<2>(0,0), Point<2>(1,1)); + + FE_Q_DG0 fe(1); + DoFHandler dofh(tr); + dofh.distribute_dofs (fe); + + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + deallog << "Total dofs=" << dofh.n_dofs() << std::endl; + + // extract constant modes and print + if (myid == 0) + { + std::vector mask(fe.n_components(), true); + + std::vector > constant_modes; + DoFTools::extract_constant_modes (dofh, mask, constant_modes); + + for (unsigned int i=0; i(); + deallog.pop(); + } + else + { + deallog.push("2d"); + test<2>(); + deallog.pop(); + } +} diff --git a/tests/mpi/extract_constant_modes_02.mpirun=4.output b/tests/mpi/extract_constant_modes_02.mpirun=4.output new file mode 100644 index 0000000000..63cce8e31f --- /dev/null +++ b/tests/mpi/extract_constant_modes_02.mpirun=4.output @@ -0,0 +1,4 @@ + +DEAL:0:2d::Total dofs=14 +DEAL:0:2d::1 1 1 1 0 +DEAL:0:2d::0 0 0 0 1 -- 2.39.5