]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add FiniteElement::get_sub_fe 5205/head
authorTimo Heister <timo.heister@gmail.com>
Fri, 6 Oct 2017 18:46:06 +0000 (14:46 -0400)
committerTimo Heister <timo.heister@gmail.com>
Fri, 20 Oct 2017 18:45:05 +0000 (14:45 -0400)
doc/news/changes/minor/20171007Heister [new file with mode: 0644]
include/deal.II/fe/fe.h
include/deal.II/fe/fe_system.h
source/fe/fe.cc
source/fe/fe_system.cc
tests/fe/get_sub_fe_01.cc [new file with mode: 0644]
tests/fe/get_sub_fe_01.debug.output [new file with mode: 0644]
tests/fe/get_sub_fe_02.cc [new file with mode: 0644]
tests/fe/get_sub_fe_02.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20171007Heister b/doc/news/changes/minor/20171007Heister
new file mode 100644 (file)
index 0000000..fc8b210
--- /dev/null
@@ -0,0 +1,3 @@
+New: FiniteElement::get_sub_fe() allows you to return a contained FiniteElement based on a ComponentMask.
+<br>
+(Timo Heister, 2017/10/07)
index e4ae7e555f0310121f7ef829c2a7fec7b810af7c..cbe84c9faaabd894328086fcb360116326ba4430 100644 (file)
@@ -1573,6 +1573,92 @@ public:
   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
@@ -2989,6 +3075,7 @@ FiniteElement<dim,spacedim>::face_system_to_component_index (const unsigned int
 
 
 
+
 template <int dim, int spacedim>
 inline
 std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
index f9218330dfdd35a2c182b011a6a228268c76ad5c..9801df566c1fd074b9f093e1a886f8d364ae969e 100644 (file)
@@ -519,6 +519,17 @@ public:
   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
index 27538dbe654f5d1215efed5ff1bb066d661c8e6e..4a5beb7e76e5b02d57b73f3c16ed1572b4dfea7d 100644 (file)
@@ -1079,6 +1079,57 @@ FiniteElement<dim,spacedim>::has_support_on_face (
 
 
 
+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
index c3b5785de0fa5e421bbf87bced55d0a597d97bb8..22e03f11f7a7fbe8a2bc0a1a408bfca9905462a5 100644 (file)
@@ -71,7 +71,6 @@ FESystem<dim,spacedim>::InternalData::~InternalData()
 }
 
 
-
 template <int dim, int spacedim>
 typename FiniteElement<dim,spacedim>::InternalDataBase &
 FESystem<dim,spacedim>::
@@ -333,6 +332,31 @@ FESystem<dim,spacedim>::clone() const
 
 
 
+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,
diff --git a/tests/fe/get_sub_fe_01.cc b/tests/fe/get_sub_fe_01.cc
new file mode 100644 (file)
index 0000000..8fedc69
--- /dev/null
@@ -0,0 +1,180 @@
+// ---------------------------------------------------------------------
+//
+// 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> ();
+}
+
diff --git a/tests/fe/get_sub_fe_01.debug.output b/tests/fe/get_sub_fe_01.debug.output
new file mode 100644 (file)
index 0000000..72d8177
--- /dev/null
@@ -0,0 +1,53 @@
+
+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)
diff --git a/tests/fe/get_sub_fe_02.cc b/tests/fe/get_sub_fe_02.cc
new file mode 100644 (file)
index 0000000..051f3ae
--- /dev/null
@@ -0,0 +1,81 @@
+// ---------------------------------------------------------------------
+//
+// 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> ();
+}
+
diff --git a/tests/fe/get_sub_fe_02.output b/tests/fe/get_sub_fe_02.output
new file mode 100644 (file)
index 0000000..b0fb89c
--- /dev/null
@@ -0,0 +1,49 @@
+
+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>

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.