--- /dev/null
+New: FiniteElement::get_sub_fe() allows you to return a contained FiniteElement based on a ComponentMask.
+<br>
+(Timo Heister, 2017/10/07)
unsigned int
element_multiplicity (const unsigned int index) const;
+ /**
+ * Return a reference to a contained finite element that matches the components
+ * selected by the given ComponentMask @p mask.
+ *
+ * For an arbitrarily nested FESystem, this function returns the inner-most
+ * FiniteElement that matches the given mask. The method fails if the @p mask
+ * does not exactly match one of the contained finite elements. It is most
+ * useful if the current object is an FESystem, as the return value can
+ * only be @p this in all other cases.
+ *
+ * Note that the returned object can be an FESystem if the
+ * mask matches it but not any of the contained objects.
+ *
+ * Let us illustrate the function with the an FESystem @p fe with 7 components:
+ * @code
+ * FESystem<2> fe_velocity(FE_Q<2>(2), 2);
+ * FE_Q<2> fe_pressure(1);
+ * FE_DGP<2> fe_dg(0);
+ * FE_BDM<2> fe_nonprim(1);
+ * FESystem<2> fe(fe_velocity, 1, fe_pressure, 1, fe_dg, 2, fe_nonprim, 1);
+ * @endcode
+ *
+ * The following table lists all possible component masks you can use:
+ * <table>
+ * <tr>
+ * <th>ComponentMask</th>
+ * <th>Result</th>
+ * <th>Description</th>
+ * </tr>
+ * <tr>
+ * <td><code>[true,true,true,true,true,true,true]</code></td>
+ * <td><code>FESystem<2>[FESystem<2>[FE_Q<2>(2)^2]-FE_Q<2>(1)-FE_DGP<2>(0)^2-FE_BDM<2>(1)]</code></td>
+ * <td>@p fe itself, the whole @p FESystem</td>
+ * </tr>
+ * <tr>
+ * <td><code>[true,true,false,false,false,false,false]</code></td>
+ * <td><code>FESystem<2>[FE_Q<2>(2)^2]</code></td>
+ * <td>just the @p fe_velocity</td>
+ * </tr>
+ * <tr>
+ * <td><code>[true,false,false,false,false,false,false]</code></td>
+ * <td><code>FE_Q<2>(2)</code></td>
+ * <td>The first component in @p fe_velocity</td>
+ * </tr>
+ * <tr>
+ * <td><code>[false,true,false,false,false,false,false]</code></td>
+ * <td><code>FE_Q<2>(2)</code></td>
+ * <td>The second component in @p fe_velocity</td>
+ * </tr>
+ * <tr>
+ * <td><code>[false,false,true,false,false,false,false]</code></td>
+ * <td><code>FE_Q<2>(1)</code></td>
+ * <td>@p fe_pressure</td>
+ * </tr>
+ * <tr>
+ * <td><code>[false,false,false,true,false,false,false]</code></td>
+ * <td><code>FE_DGP<2>(0)</code></td>
+ * <td>first copy of @p fe_dg</td>
+ * </tr>
+ * <tr>
+ * <td><code>[false,false,false,false,true,false,false]</code></td>
+ * <td><code>FE_DGP<2>(0)</code></td>
+ * <td>second copy of @p fe_dg</td>
+ * </tr>
+ * <tr>
+ * <td><code>[false,false,false,false,false,true,true]</code></td>
+ * <td><code>FE_BDM<2>(1)</code></td>
+ * <td>both components of @p fe_nonprim</td>
+ * </tr>
+ * </table>
+ */
+ const FiniteElement<dim,spacedim> &
+ get_sub_fe (const ComponentMask &mask) const;
+
+ /**
+ * Return a reference to a contained finite element that matches the components
+ * @p n_selected_components components starting at component with index
+ * @p first_component.
+ *
+ * See the other get_sub_fe() function above for more details.
+ */
+ virtual
+ const FiniteElement<dim,spacedim> &
+ get_sub_fe (const unsigned int first_component,
+ const unsigned int n_selected_components) const;
+
/**
* Return for shape function @p index the base element it belongs to, the
* number of the copy of this base element (which is between zero and the
+
template <int dim, int spacedim>
inline
std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
UpdateFlags
requires_update_flags (const UpdateFlags update_flags) const;
+ // make variant with ComponentMask also available:
+ using FiniteElement<dim,spacedim>::get_sub_fe;
+
+ /**
+ * @copydoc FiniteElement<dim,spacedim>::get_sub_fe()
+ */
+ virtual
+ const FiniteElement<dim,spacedim> &
+ get_sub_fe (const unsigned int first_component,
+ const unsigned int n_selected_components) const;
+
/**
* Return the value of the @p ith shape function at the point @p p. @p p is
* a point on the reference element. Since this finite element is always
+template <int dim, int spacedim>
+const FiniteElement<dim,spacedim> &
+FiniteElement<dim,spacedim>::
+get_sub_fe (const ComponentMask &mask) const
+{
+ // Translate the ComponentMask into first_selected and n_components after
+ // some error checking:
+ const unsigned int n_total_components = this->n_components();
+ Assert((n_total_components == mask.size()) || (mask.size()==0),
+ ExcMessage("The given ComponentMask has the wrong size."));
+
+ const unsigned int n_selected = mask.n_selected_components(n_total_components);
+ Assert(n_selected>0,
+ ExcMessage("You need at least one selected component."));
+
+ const unsigned int first_selected = mask.first_selected_component(n_total_components);
+
+#ifdef DEBUG
+ // check that it is contiguous:
+ for (unsigned int c=0; c<n_total_components; ++c)
+ Assert((c<first_selected && (!mask[c]))
+ ||
+ (c>=first_selected && c<first_selected+n_selected && mask[c])
+ ||
+ (c>=first_selected+n_selected && !mask[c]),
+ ExcMessage("Error: the given ComponentMask is not contiguous!"));
+#endif
+
+ return get_sub_fe(first_selected, n_selected);
+}
+
+
+
+template <int dim, int spacedim>
+const FiniteElement<dim,spacedim> &
+FiniteElement<dim,spacedim>::
+get_sub_fe (const unsigned int first_component,
+ const unsigned int n_selected_components) const
+{
+ // No complicated logic is needed here, because it is overridden in
+ // FESystem<dim,spacedim>. Just make sure that what the user chose is valid:
+ Assert(first_component == 0 && n_selected_components == this->n_components(),
+ ExcMessage("You can only select a whole FiniteElement, not a part of one."));
+
+ (void)first_component;
+ (void)n_selected_components;
+ return *this;
+}
+
+
+
template <int dim, int spacedim>
std::pair<Table<2,bool>, std::vector<unsigned int> >
FiniteElement<dim,spacedim>::get_constant_modes () const
}
-
template <int dim, int spacedim>
typename FiniteElement<dim,spacedim>::InternalDataBase &
FESystem<dim,spacedim>::
+template <int dim, int spacedim>
+const FiniteElement<dim,spacedim> &
+FESystem<dim,spacedim>::
+get_sub_fe (const unsigned int first_component,
+ const unsigned int n_selected_components) const
+{
+ Assert(first_component+n_selected_components <= this->n_components(),
+ ExcMessage("Invalid arguments (not a part of this FiniteElement)."));
+
+ const unsigned int base_index = this->component_to_base_table[first_component].first.first;
+ const unsigned int component_in_base = this->component_to_base_table[first_component].first.second;
+ const unsigned int base_components = this->base_element(base_index).n_components();
+
+ // Only select our child base_index if that is all the user wanted. Error
+ // handling will be done inside the recursion.
+ if (n_selected_components<=base_components)
+ return this->base_element(base_index).get_sub_fe(component_in_base, n_selected_components);
+
+ Assert(n_selected_components == this->n_components(),
+ ExcMessage("You can not select a part of a FiniteElement."));
+ return *this;
+}
+
+
+
template <int dim, int spacedim>
double
FESystem<dim,spacedim>::shape_value (const unsigned int i,
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2002 - 2015 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 FiniteElement::get_sub_fe()
+
+#include "../tests.h"
+#include <iostream>
+
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_bdm.h>
+#include <deal.II/fe/fe_system.h>
+
+template <int dim>
+void works(const FiniteElement<dim> &fe, const ComponentMask &m)
+{
+ deallog << "FE: " << fe.get_name()
+ << " mask: " << m << std::endl;
+
+ const FiniteElement<dim> &child = fe.get_sub_fe(m);
+
+ deallog << " worked: " << child.get_name() << std::endl;
+}
+
+template <int dim>
+void fails(const FiniteElement<dim> &fe, const ComponentMask &m)
+{
+ deallog << "FE: " << fe.get_name()
+ << " mask: " << m << std::endl;
+
+ try
+ {
+ const FiniteElement<dim> &child = fe.get_sub_fe(m);
+ deallog << " ERROR: we succeeded and got " << child.get_name()
+ << " but we should have failed!" << std::endl;
+ }
+ catch (...)
+ {
+ deallog << " failed as expected" << std::endl;
+ }
+}
+
+
+template <int dim>
+void check ()
+{
+ auto mask_none = [] (const unsigned int n_components) -> ComponentMask
+ {
+ return ComponentMask(n_components, false);
+ };
+ auto mask_all = [] (const unsigned int n_components) -> ComponentMask
+ {
+ return ComponentMask(n_components, true);
+ };
+ auto mask_single = [] (const unsigned int n_components,
+ const unsigned int single) -> ComponentMask
+ {
+ ComponentMask c(n_components, false);
+ c.set(single, true);
+ return c;
+ };
+ auto mask = [] (const unsigned int n_components,
+ const unsigned int first,
+ const unsigned int last) -> ComponentMask
+ {
+ ComponentMask c(n_components, false);
+ for (unsigned int i=first; i<=last; ++i)
+ c.set(i, true);
+ return c;
+ };
+
+ FE_Q<dim> fe_q(2);
+ FE_DGP<dim> fe_dg(2);
+ FE_BDM<dim> fe_nonprim(1);
+
+ // simple FE:
+ {
+ works(fe_q, mask_all(1));
+ fails(fe_q, mask_none(1));
+ }
+
+ // simple non-primitive:
+ {
+ works(fe_nonprim, ComponentMask()); // un-sized "select all"
+ works(fe_nonprim, mask_all(dim));
+ fails(fe_nonprim, mask_none(dim));
+ fails(fe_nonprim, mask_single(dim, 0)); // can not select part of it
+ }
+
+ // remove system:
+ {
+ FESystem<dim> fe(fe_q, 1);
+ FESystem<dim> fe2(fe, 1);
+ works(fe, mask_all(1));
+ works(fe2, mask_all(1));
+ }
+
+ // part of system:
+ {
+ FESystem<dim> fe(fe_q, 3);
+ works(fe, mask_all(3)); // keep system
+ works(fe, mask_single(3, 1)); // select one component
+ fails(fe, mask(3, 1, 2)); // select more than one component but not the whole FESystem
+ }
+
+ // more systems:
+ {
+ FESystem<dim> fe(fe_nonprim,1,fe_dg,1);
+ works(fe, mask(dim+1, 0, dim)); // nonprimitive
+ works(fe, mask_single(dim+1, dim)); // select one component
+ }
+ {
+ FESystem<dim> fe(fe_nonprim,1,fe_dg,2);
+ // non-contiguous is a fail!
+ auto m = mask(dim+2, 0, dim);
+ m.set(1, false);
+ m.set(dim+1,true);
+ fails(fe, m);
+ }
+ {
+ FESystem<dim> fe(fe_q,dim,fe_q,1);
+ works(fe, mask_single(dim+1, dim));
+ fails(fe, mask(dim+1, 0, dim-1)); // can not select the dim fe_q
+ }
+ {
+ FESystem<dim> qs(fe_q,dim);
+ FESystem<dim> fe(qs,1,fe_q,1);
+ works(fe, mask(dim+1, 0, dim-1)); // but works if nested
+ }
+
+ {
+ FESystem<dim> fe(fe_q,2,fe_nonprim,2,fe_dg,1);
+ works(fe, mask_single(2*dim+3, 0));
+ works(fe, mask_single(2*dim+3, 1));
+ works(fe, mask_single(2*dim+3, 2*dim+3-1));
+ works(fe, mask(2*dim+3, 2, 2+dim-1)); // 1st nonprimitive
+ works(fe, mask(2*dim+3, 2+dim, 2+2*dim-1)); // 2nd nonprimitive
+ }
+
+ // nesting doll:
+ {
+ FESystem<dim> inner(fe_dg,1);
+ FESystem<dim> fe(fe_nonprim,1,inner,2);
+ FESystem<dim> outer(fe,2);
+ FESystem<dim> outer_outer(outer,1);
+
+ works(fe, mask_single(dim+2, dim));
+ works(fe, mask_single(dim+2, dim+1));
+ works(outer, mask_single(2*(dim+2), dim));
+ works(outer_outer, mask_single(2*(dim+2), dim));
+ }
+
+}
+
+
+
+
+int main ()
+{
+ initlog();
+ deal_II_exceptions::disable_abort_on_exception();
+
+ check<2> ();
+}
+
--- /dev/null
+
+DEAL::FE: FE_Q<2>(2) mask: [true]
+DEAL:: worked: FE_Q<2>(2)
+DEAL::FE: FE_Q<2>(2) mask: [false]
+DEAL:: failed as expected
+DEAL::FE: FE_BDM<2>(1) mask: [all components selected]
+DEAL:: worked: FE_BDM<2>(1)
+DEAL::FE: FE_BDM<2>(1) mask: [true,true]
+DEAL:: worked: FE_BDM<2>(1)
+DEAL::FE: FE_BDM<2>(1) mask: [false,false]
+DEAL:: failed as expected
+DEAL::FE: FE_BDM<2>(1) mask: [true,false]
+DEAL:: failed as expected
+DEAL::FE: FESystem<2>[FE_Q<2>(2)] mask: [true]
+DEAL:: worked: FE_Q<2>(2)
+DEAL::FE: FESystem<2>[FESystem<2>[FE_Q<2>(2)]] mask: [true]
+DEAL:: worked: FE_Q<2>(2)
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^3] mask: [true,true,true]
+DEAL:: worked: FESystem<2>[FE_Q<2>(2)^3]
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^3] mask: [false,true,false]
+DEAL:: worked: FE_Q<2>(2)
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^3] mask: [false,true,true]
+DEAL:: failed as expected
+DEAL::FE: FESystem<2>[FE_BDM<2>(1)-FE_DGP<2>(2)] mask: [true,true,true]
+DEAL:: worked: FESystem<2>[FE_BDM<2>(1)-FE_DGP<2>(2)]
+DEAL::FE: FESystem<2>[FE_BDM<2>(1)-FE_DGP<2>(2)] mask: [false,false,true]
+DEAL:: worked: FE_DGP<2>(2)
+DEAL::FE: FESystem<2>[FE_BDM<2>(1)-FE_DGP<2>(2)^2] mask: [true,false,true,true]
+DEAL:: failed as expected
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(2)] mask: [false,false,true]
+DEAL:: worked: FE_Q<2>(2)
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(2)] mask: [true,true,false]
+DEAL:: failed as expected
+DEAL::FE: FESystem<2>[FESystem<2>[FE_Q<2>(2)^2]-FE_Q<2>(2)] mask: [true,true,false]
+DEAL:: worked: FESystem<2>[FE_Q<2>(2)^2]
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^2-FE_BDM<2>(1)^2-FE_DGP<2>(2)] mask: [true,false,false,false,false,false,false]
+DEAL:: worked: FE_Q<2>(2)
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^2-FE_BDM<2>(1)^2-FE_DGP<2>(2)] mask: [false,true,false,false,false,false,false]
+DEAL:: worked: FE_Q<2>(2)
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^2-FE_BDM<2>(1)^2-FE_DGP<2>(2)] mask: [false,false,false,false,false,false,true]
+DEAL:: worked: FE_DGP<2>(2)
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^2-FE_BDM<2>(1)^2-FE_DGP<2>(2)] mask: [false,false,true,true,false,false,false]
+DEAL:: worked: FE_BDM<2>(1)
+DEAL::FE: FESystem<2>[FE_Q<2>(2)^2-FE_BDM<2>(1)^2-FE_DGP<2>(2)] mask: [false,false,false,false,true,true,false]
+DEAL:: worked: FE_BDM<2>(1)
+DEAL::FE: FESystem<2>[FE_BDM<2>(1)-FESystem<2>[FE_DGP<2>(2)]^2] mask: [false,false,true,false]
+DEAL:: worked: FE_DGP<2>(2)
+DEAL::FE: FESystem<2>[FE_BDM<2>(1)-FESystem<2>[FE_DGP<2>(2)]^2] mask: [false,false,false,true]
+DEAL:: worked: FE_DGP<2>(2)
+DEAL::FE: FESystem<2>[FESystem<2>[FE_BDM<2>(1)-FESystem<2>[FE_DGP<2>(2)]^2]^2] mask: [false,false,true,false,false,false,false,false]
+DEAL:: worked: FE_DGP<2>(2)
+DEAL::FE: FESystem<2>[FESystem<2>[FESystem<2>[FE_BDM<2>(1)-FESystem<2>[FE_DGP<2>(2)]^2]^2]] mask: [false,false,true,false,false,false,false,false]
+DEAL:: worked: FE_DGP<2>(2)
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2002 - 2015 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 FiniteElement::get_sub_fe(), example used in documentation
+
+#include "../tests.h"
+#include <iostream>
+
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_bdm.h>
+#include <deal.II/fe/fe_system.h>
+
+template <int dim>
+void check ()
+{
+ FESystem<2> fe_velocity(FE_Q<2>(2), 2);
+ FE_Q<2> fe_pressure(1);
+ FE_DGP<2> fe_dg(0);
+ FE_BDM<2> fe_nonprim(1);
+
+ FESystem<2> fe(fe_velocity, 1, fe_pressure, 1, fe_dg, 2, fe_nonprim, 1);
+
+ // same using component masks to copy over:
+ auto run = [&](const unsigned int first,
+ const unsigned int n,
+ const std::string &desc)
+ {
+ const unsigned int n_components = fe.n_components();
+
+ ComponentMask mask(n_components, false);
+ for (unsigned int i=first; i<first+n; ++i)
+ mask.set(i, true);
+
+ deallog
+ << "<tr>\n"
+ << "<td><code>" << mask << "</code></td>\n"
+ << "<td><code>" << fe.get_sub_fe(mask).get_name() << "</code></td>\n"
+ << "<td>" << desc << "</td>\n"
+ << "</tr>\n";
+
+ // we should be able to use ComponentMask or indices:
+ Assert(fe.get_sub_fe(mask) == fe.get_sub_fe(first, n),
+ ExcInternalError());
+ };
+
+ deallog << "\n<table>\n"
+ << "<tr>\n<th>ComponentMask</th>\n<th>Result</th>\n<th>Description</th>\n</tr>\n";
+ run(0, 7, "@p fe itself, the whole @p FESystem");
+ run(0, 2, "just the @p fe_velocity");
+ run(0, 1, "The first component in @p fe_velocity");
+ run(1, 1, "The second component in @p fe_velocity");
+ run(2, 1, "@p fe_pressure");
+ run(3, 1, "first copy of @p fe_dg");
+ run(4, 1, "second copy of @p fe_dg");
+ run(5, 2, "both components of @p fe_nonprim");
+ deallog << "</table>" << std::endl;
+}
+
+int main ()
+{
+ initlog();
+
+ check<2> ();
+}
+
--- /dev/null
+
+DEAL::
+<table>
+<tr>
+<th>ComponentMask</th>
+<th>Result</th>
+<th>Description</th>
+</tr>
+<tr>
+<td><code>[true,true,true,true,true,true,true]</code></td>
+<td><code>FESystem<2>[FESystem<2>[FE_Q<2>(2)^2]-FE_Q<2>(1)-FE_DGP<2>(0)^2-FE_BDM<2>(1)]</code></td>
+<td>@p fe itself, the whole @p FESystem</td>
+</tr>
+<tr>
+<td><code>[true,true,false,false,false,false,false]</code></td>
+<td><code>FESystem<2>[FE_Q<2>(2)^2]</code></td>
+<td>just the @p fe_velocity</td>
+</tr>
+<tr>
+<td><code>[true,false,false,false,false,false,false]</code></td>
+<td><code>FE_Q<2>(2)</code></td>
+<td>The first component in @p fe_velocity</td>
+</tr>
+<tr>
+<td><code>[false,true,false,false,false,false,false]</code></td>
+<td><code>FE_Q<2>(2)</code></td>
+<td>The second component in @p fe_velocity</td>
+</tr>
+<tr>
+<td><code>[false,false,true,false,false,false,false]</code></td>
+<td><code>FE_Q<2>(1)</code></td>
+<td>@p fe_pressure</td>
+</tr>
+<tr>
+<td><code>[false,false,false,true,false,false,false]</code></td>
+<td><code>FE_DGP<2>(0)</code></td>
+<td>first copy of @p fe_dg</td>
+</tr>
+<tr>
+<td><code>[false,false,false,false,true,false,false]</code></td>
+<td><code>FE_DGP<2>(0)</code></td>
+<td>second copy of @p fe_dg</td>
+</tr>
+<tr>
+<td><code>[false,false,false,false,false,true,true]</code></td>
+<td><code>FE_BDM<2>(1)</code></td>
+<td>both components of @p fe_nonprim</td>
+</tr>
+</table>
projects / dealii.git / commitdiff