<h3>Specific improvements</h3>
<ol>
+ <li> 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.
+ <br>
+ (Martin Kronbichler, 2014/04/04)
+ </li>
+
<li> 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:
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<Table<2,bool>,std::vector<unsigned int> >
+ get_constant_modes () const;
//@}
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
protected:
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
/**
* Determine an estimate for the
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
private:
/**
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
protected:
virtual
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
private:
/**
/**
* Returns a list of constant modes of the element.
*/
- virtual Table<2,bool> get_constant_modes () const;
+ virtual std::pair<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
virtual std::size_t memory_consumption () const;
virtual FiniteElement<dim> *clone() const;
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
/**
* @name Functions to support hp
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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
+
protected:
/**
* @p clone function instead of a copy constructor.
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
protected:
/**
const VectorSlice<const std::vector<std::vector<double> > > &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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
virtual std::size_t memory_consumption () const;
virtual FiniteElement<dim> *clone() const;
* 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<Table<2,bool>, std::vector<unsigned int> >
+ get_constant_modes () const;
/**
* @name Functions to support hp
Assert (component_mask.represents_n_components(n_components),
ExcDimensionMismatch(n_components,
component_mask.size()));
- std::vector<unsigned int> localized_component (n_components,
- numbers::invalid_unsigned_int);
- unsigned int n_components_selected = 0;
- for (unsigned int i=0; i<n_components; ++i)
- if (component_mask[i] == true)
- localized_component[i] = n_components_selected++;
std::vector<unsigned char> dofs_by_component (dof_handler.n_locally_owned_dofs());
internal::get_component_association (dof_handler, component_mask,
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<bool>(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<DH::dimension,DH::space_dimension>
fe_collection (dof_handler.get_fe());
std::vector<Table<2,bool> > element_constant_modes;
+ std::vector<std::vector<std::pair<unsigned int, unsigned int> > >
+ constant_mode_to_component_translation(n_components);
+ unsigned int n_constant_modes = 0;
for (unsigned int f=0; f<fe_collection.size(); ++f)
- element_constant_modes.push_back(fe_collection[f].get_constant_modes());
+ {
+ std::pair<Table<2,bool>, std::vector<unsigned int> > data
+ = fe_collection[f].get_constant_modes();
+ element_constant_modes.push_back(data.first);
+ if (f==0)
+ for (unsigned int i=0; i<data.second.size(); ++i)
+ if (component_mask[data.second[i]])
+ constant_mode_to_component_translation[data.second[i]].
+ push_back(std::make_pair(n_constant_modes++,i));
+ AssertDimension(element_constant_modes.back().n_rows(),
+ element_constant_modes[0].n_rows());
+ }
+
+ // First count the number of dofs in the current component.
+ constant_modes.clear ();
+ constant_modes.resize (n_constant_modes, std::vector<bool>(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();
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; j<constant_mode_to_component_translation[comp].size(); ++j)
+ constant_modes[constant_mode_to_component_translation[comp][j].first]
+ [component_numbering[loc_index]] =
+ element_constant_modes[cell->active_fe_index()]
+ (constant_mode_to_component_translation[comp][j].second,i);
}
}
}
#include <deal.II/fe/fe_dgq.h>
#include <deal.II/fe/fe_dgp.h>
#include <deal.II/fe/fe_face.h>
+#include <deal.II/fe/fe_q_dg0.h>
#include <deal.II/fe/fe_q_hierarchical.h>
#include <deal.II/fe/fe_system.h>
#include <fstream>
template<int dim>
void print_constant_modes(const FiniteElement<dim> &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<constant_modes.n_rows(); ++r)
{
for (unsigned int c=0; c<constant_modes.n_cols(); ++c)
print_constant_modes(FE_FaceP<dim>(1));
print_constant_modes(FESystem<dim>(FE_Q<dim>(1), 2, FE_Q<dim>(2), 1));
print_constant_modes(FESystem<dim>(FE_DGP<dim>(1), 1, FE_Q_iso_Q1<dim>(2), 1));
+ print_constant_modes(FE_Q_DG0<dim>(1));
+ print_constant_modes(FESystem<dim>(FE_Q_DG0<dim>(2), 1, FE_Q<dim>(1), 2));
+ print_constant_modes(FESystem<dim>(FE_Q<dim>(1), 2, FE_Q_DG0<dim>(1), 2));
}
template <>
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::
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::
--- /dev/null
+// ---------------------------------------------------------------------
+// $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 <deal.II/base/logstream.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/dofs/dof_renumbering.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/dofs/dof_accessor.h>
+
+#include <deal.II/fe/fe_q_dg0.h>
+
+#include <fstream>
+
+
+template<int dim>
+void test()
+{
+ unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+ parallel::distributed::Triangulation<dim> tr(MPI_COMM_WORLD);
+
+ std::vector<unsigned int> sub(2);
+ sub[0] = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);
+ sub[1] = 1;
+ GridGenerator::subdivided_hyper_rectangle(static_cast<Triangulation<dim>&>(tr),
+ sub, Point<2>(0,0), Point<2>(1,1));
+
+ FE_Q_DG0<dim> fe(1);
+ DoFHandler<dim> 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<bool> mask(fe.n_components(), true);
+
+ std::vector<std::vector<bool> > constant_modes;
+ DoFTools::extract_constant_modes (dofh, mask, constant_modes);
+
+ for (unsigned int i=0; i<constant_modes.size(); ++i)
+ {
+ for (unsigned int j=0; j<constant_modes[i].size(); ++j)
+ deallog << (constant_modes[i][j] ? '1' : '0') << ' ';
+ deallog << std::endl;
+ }
+ }
+}
+
+
+int main(int argc, char *argv[])
+{
+ Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+ unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+
+ deallog.push(Utilities::int_to_string(myid));
+
+ if (myid == 0)
+ {
+ std::ofstream logfile("output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ deallog.push("2d");
+ test<2>();
+ deallog.pop();
+ }
+ else
+ {
+ deallog.push("2d");
+ test<2>();
+ deallog.pop();
+ }
+}
--- /dev/null
+
+DEAL:0:2d::Total dofs=14
+DEAL:0:2d::1 1 1 1 0
+DEAL:0:2d::0 0 0 0 1