From: Denis Davydov Date: Tue, 29 Aug 2017 19:35:55 +0000 (+0200) Subject: explain FiniteElement::system_to_xyz by considering Q2xQ2xQ1 example X-Git-Tag: v9.0.0-rc1~1141^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a5f4df4b54733edb36cf676c6141d185606fbaa3;p=dealii.git explain FiniteElement::system_to_xyz by considering Q2xQ2xQ1 example --- diff --git a/doc/doxygen/images/fe_system_example.png b/doc/doxygen/images/fe_system_example.png new file mode 100644 index 0000000000..467ddb0d9c Binary files /dev/null and b/doc/doxygen/images/fe_system_example.png differ diff --git a/include/deal.II/fe/fe.h b/include/deal.II/fe/fe.h index 165d794211..ebf915637c 100644 --- a/include/deal.II/fe/fe.h +++ b/include/deal.II/fe/fe.h @@ -149,6 +149,67 @@ template 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 fe_basis(FE_Q(2), dim, FE_Q(1),1); + * @endcode + * with dim=2. 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. * *

Support points

* @@ -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 index 0000000000..52376900a6 --- /dev/null +++ b/tests/fe/fe_system_components_base.cc @@ -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 +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +template +void test(const bool renumber = false) +{ + Triangulation triangulation; + FESystem fe_basis(FE_Q(2), dim, FE_Q(1),1); + DoFHandler dof_handler (triangulation); + GridGenerator::hyper_cube(triangulation); + dof_handler.distribute_dofs(fe_basis); + + if (renumber) + { + std::vector 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 > support_points; + MappingQ1 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 index 0000000000..dda6e8d8ad --- /dev/null +++ b/tests/fe/fe_system_components_base.output @@ -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 |