]> https://gitweb.dealii.org/ - dealii.git/commitdiff
explain FiniteElement::system_to_xyz by considering Q2xQ2xQ1 example 4990/head
authorDenis Davydov <davydden@gmail.com>
Tue, 29 Aug 2017 19:35:55 +0000 (21:35 +0200)
committerDenis Davydov <davydden@gmail.com>
Thu, 31 Aug 2017 07:17:09 +0000 (09:17 +0200)
doc/doxygen/images/fe_system_example.png [new file with mode: 0644]
include/deal.II/fe/fe.h
tests/fe/fe_system_components_base.cc [new file with mode: 0644]
tests/fe/fe_system_components_base.output [new file with mode: 0644]

diff --git a/doc/doxygen/images/fe_system_example.png b/doc/doxygen/images/fe_system_example.png
new file mode 100644 (file)
index 0000000..467ddb0
Binary files /dev/null and b/doc/doxygen/images/fe_system_example.png differ
index 165d794211d5f3d6df762ad7babc9dace5dd4335..ebf915637c165b03d898cd2194f7381092befd23 100644 (file)
@@ -149,6 +149,67 @@ template <int dim, int spacedim> class FESystem;
  * FiniteElement::system_to_block_index(). The number of blocks of a finite
  * element can be determined by FiniteElement::n_blocks().
  *
+ * To better illustrate these concepts, let's consider the following
+ * example of the multi-component system
+ * @code
+ * FESystem<dim> fe_basis(FE_Q<dim>(2), dim, FE_Q<dim>(1),1);
+ * @endcode
+ * with <code>dim=2</code>. The resulting finite element has 3 components:
+ * two that come from the quadratic element and one from the linear element.
+ * If, for example, this system were used to discretize a problem in fluid
+ * dynamics then one could think of the first two components representing a
+ * vector-valued velocity field whereas the last one corresponds to the scalar
+ * pressure field. Without degree-of-freedom (DoF) renumbering this finite element will
+ * produce the following distribution of local DoFs:
+ *
+ * @image html fe_system_example.png DoF indices
+ *
+ * Using the two functions FiniteElement::system_to_component_index() and
+ * FiniteElement::system_to_base_index() one can get the
+ * following information for each degree-of-freedom "i":
+ * @code
+ * const unsigned int component     = fe_basis.system_to_component_index(i).first;
+ * const unsigned int within_base   = fe_basis.system_to_component_index(i).second;
+ * const unsigned int base          = fe_basis.system_to_base_index(i).first.first;
+ * const unsigned int multiplicity  = fe_basis.system_to_base_index(i).first.second;
+ * const unsigned int within_base_  = fe_basis.system_to_base_index(i).second; // same as above
+ * @endcode
+ * which will result in:
+ *
+ * | DoF    | Component  | Base element | Shape function within base | Multiplicity |
+ * | :----: | :--------: | :----------: | :------------------------: | :----------: |
+ * |      0 |          0 |            0 |                          0 |            0 |
+ * |      1 |          1 |            0 |                          0 |            1 |
+ * |      2 |          2 |            1 |                          0 |            0 |
+ * |      3 |          0 |            0 |                          1 |            0 |
+ * |      4 |          1 |            0 |                          1 |            1 |
+ * |      5 |          2 |            1 |                          1 |            0 |
+ * |      6 |          0 |            0 |                          2 |            0 |
+ * |      7 |          1 |            0 |                          2 |            1 |
+ * |      8 |          2 |            1 |                          2 |            0 |
+ * |      9 |          0 |            0 |                          3 |            0 |
+ * |     10 |          1 |            0 |                          3 |            1 |
+ * |     11 |          2 |            1 |                          3 |            0 |
+ * |     12 |          0 |            0 |                          4 |            0 |
+ * |     13 |          1 |            0 |                          4 |            1 |
+ * |     14 |          0 |            0 |                          5 |            0 |
+ * |     15 |          1 |            0 |                          5 |            1 |
+ * |     16 |          0 |            0 |                          6 |            0 |
+ * |     17 |          1 |            0 |                          6 |            1 |
+ * |     18 |          0 |            0 |                          7 |            0 |
+ * |     19 |          1 |            0 |                          7 |            1 |
+ * |     20 |          0 |            0 |                          8 |            0 |
+ * |     21 |          1 |            0 |                          8 |            1 |
+ *
+ * What we see is the following: there are a total of 22 degrees-of-freedom on this
+ * element with components ranging from 0 to 2. Each DoF corresponds to
+ * one of the two base elements used to build FESystem : $\mathbb Q_2$ or $\mathbb Q_1$.
+ * Since FE_Q are primitive elements, we have a total of 9 distinct
+ * scalar-valued shape functions for the quadratic element and 4 for the linear element.
+ * Finally, for DoFs corresponding to the first base element multiplicity
+ * is either zero or one, meaning that we use the same scalar valued $\mathbb Q_2$
+ * for both $x$ and $y$ components of the velocity field $\mathbb Q_2 \otimes \mathbb Q_2$.
+ * For DoFs corresponding to the second base element multiplicity is zero.
  *
  * <h4>Support points</h4>
  *
@@ -1287,6 +1348,8 @@ public:
    * than one vector-component). For this information, refer to the
    * #system_to_base_table field and the system_to_base_index() function.
    *
+   * See the class description above for an example of how this function is typically used.
+   *
    * The use of this function is explained extensively in the step-8 and
    * @ref step_20 "step-20"
    * tutorial programs as well as in the
@@ -1504,6 +1567,8 @@ public:
    * composed of other elements and at least one of them is vector-valued
    * itself.
    *
+   * See the class documentation above for an example of how this function is typically used.
+   *
    * This function returns valid values also in the case of vector-valued
    * (i.e. non-primitive) shape functions, in contrast to the
    * system_to_component_index() function.
diff --git a/tests/fe/fe_system_components_base.cc b/tests/fe/fe_system_components_base.cc
new file mode 100644 (file)
index 0000000..5237690
--- /dev/null
@@ -0,0 +1,110 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2003 - 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.
+//
+// ---------------------------------------------------------------------
+
+// Check FESystem system_to_component_index() and system_to_base_index()
+// for a simple FE_Q case and prints the results.
+// If at some point these results will change (due to internal refactoring),
+// update the documentation of FiniteElement::system_to_base_index()
+
+#include "../tests.h"
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/fe_tools.h>
+#include <deal.II/fe/mapping_q1.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/dofs/dof_renumbering.h>
+
+#include <string>
+#include <iomanip>
+
+template <int dim>
+void test(const bool renumber = false)
+{
+  Triangulation<dim> triangulation;
+  FESystem<dim> fe_basis(FE_Q<dim>(2), dim, FE_Q<dim>(1),1);
+  DoFHandler<dim> dof_handler (triangulation);
+  GridGenerator::hyper_cube(triangulation);
+  dof_handler.distribute_dofs(fe_basis);
+
+  if (renumber)
+    {
+      std::vector<unsigned int> component_to_block_indices(dim+1);
+      for (int i = 0; i < dim; ++i)
+        component_to_block_indices[i] = 0;
+      component_to_block_indices[dim] = 1;
+      DoFRenumbering::component_wise (dof_handler, component_to_block_indices);
+    }
+
+  const unsigned int n_dofs = dof_handler.n_dofs();
+
+  std::cout << " * | DoF    | Component  | Base element | Shape function within base | Multiplicity |" << std::endl
+            << " * | :----: | :--------: | :----------: | :------------------------: | :----------: |" << std::endl;
+
+  for (unsigned int i = 0; i < n_dofs; ++i)
+    {
+      const unsigned int component     = fe_basis.system_to_component_index(i).first;
+      const unsigned int within_base   = fe_basis.system_to_component_index(i).second;
+      const unsigned int base          = fe_basis.system_to_base_index(i).first.first;
+      const unsigned int multiplicity  = fe_basis.system_to_base_index(i).first.second;
+      const unsigned int within_base_  = fe_basis.system_to_base_index(i).second; // same as above
+      std::cout << std::setfill(' ')
+                << " * | " << std::setw(6) << i
+                << " | " << std::setw(10) << component
+                << " | " << std::setw(12) << base
+                << " | " << std::setw(26) << within_base
+                << " | " << std::setw(12) << multiplicity
+                << " |" << std::endl;
+    }
+
+  // print grid and DoFs for visual inspection
+  if (true)
+    {
+      std::map<types::global_dof_index, Point<dim> > support_points;
+      MappingQ1<dim> mapping;
+      DoFTools::map_dofs_to_support_points(mapping, dof_handler, support_points);
+
+      const std::string filename =
+        "grid" + Utilities::int_to_string(dim) + Utilities::int_to_string(renumber) + ".gp";
+      std::ofstream f(filename.c_str());
+
+      f << "set terminal png size 420,440 enhanced font \"Helvetica,16\"" << std::endl
+        << "set output \"grid" << Utilities::int_to_string(dim) << Utilities::int_to_string(renumber) << ".png\"" << std::endl
+        << "set size square" << std::endl
+        << "set view equal xy" << std::endl
+        << "unset xtics" << std::endl
+        << "unset ytics" << std::endl
+        << "unset border" << std::endl
+        << "set xrange [0: 1.05]" << std::endl
+        << "set yrange [0: 1.05]" << std::endl
+        << "plot '-' using 1:2 with lines notitle, '-' with labels point pt 2 offset 0.5,0.5 notitle" << std::endl;
+      GridOut().write_gnuplot (triangulation, f);
+      f << "e" << std::endl;
+
+      DoFTools::write_gnuplot_dof_support_point_info(f,
+                                                     support_points);
+      f << "e" << std::endl;
+    }
+}
+
+int
+main()
+{
+  test<2>(false);
+
+  return 0;
+}
diff --git a/tests/fe/fe_system_components_base.output b/tests/fe/fe_system_components_base.output
new file mode 100644 (file)
index 0000000..dda6e8d
--- /dev/null
@@ -0,0 +1,24 @@
+ * | DoF    | Component  | Base element | Shape function within base | Multiplicity |
+ * | :----: | :--------: | :----------: | :------------------------: | :----------: |
+ * |      0 |          0 |            0 |                          0 |            0 |
+ * |      1 |          1 |            0 |                          0 |            1 |
+ * |      2 |          2 |            1 |                          0 |            0 |
+ * |      3 |          0 |            0 |                          1 |            0 |
+ * |      4 |          1 |            0 |                          1 |            1 |
+ * |      5 |          2 |            1 |                          1 |            0 |
+ * |      6 |          0 |            0 |                          2 |            0 |
+ * |      7 |          1 |            0 |                          2 |            1 |
+ * |      8 |          2 |            1 |                          2 |            0 |
+ * |      9 |          0 |            0 |                          3 |            0 |
+ * |     10 |          1 |            0 |                          3 |            1 |
+ * |     11 |          2 |            1 |                          3 |            0 |
+ * |     12 |          0 |            0 |                          4 |            0 |
+ * |     13 |          1 |            0 |                          4 |            1 |
+ * |     14 |          0 |            0 |                          5 |            0 |
+ * |     15 |          1 |            0 |                          5 |            1 |
+ * |     16 |          0 |            0 |                          6 |            0 |
+ * |     17 |          1 |            0 |                          6 |            1 |
+ * |     18 |          0 |            0 |                          7 |            0 |
+ * |     19 |          1 |            0 |                          7 |            1 |
+ * |     20 |          0 |            0 |                          8 |            0 |
+ * |     21 |          1 |            0 |                          8 |            1 |

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.