From a5f4df4b54733edb36cf676c6141d185606fbaa3 Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Tue, 29 Aug 2017 21:35:55 +0200 Subject: [PATCH] explain FiniteElement::system_to_xyz by considering Q2xQ2xQ1 example --- doc/doxygen/images/fe_system_example.png | Bin 0 -> 3397 bytes include/deal.II/fe/fe.h | 65 +++++++++++++ tests/fe/fe_system_components_base.cc | 110 ++++++++++++++++++++++ tests/fe/fe_system_components_base.output | 24 +++++ 4 files changed, 199 insertions(+) create mode 100644 doc/doxygen/images/fe_system_example.png create mode 100644 tests/fe/fe_system_components_base.cc create mode 100644 tests/fe/fe_system_components_base.output diff --git a/doc/doxygen/images/fe_system_example.png b/doc/doxygen/images/fe_system_example.png new file mode 100644 index 0000000000000000000000000000000000000000..467ddb0d9c42e50e55c73b7e0a3e15dd28667118 GIT binary patch literal 3397 zcmeHJdpOhkAOCiHej%wy8KPox7mi6px|wWboV8d}Zn@uLF2jx!PN^`LHOi$Hhtp7U zZ%lIuTRJ!x*_1}P%+SPKzQ5IfzvuZq=dbhU`8@CU=ktEOKF|B}zWniiKDV8a7i6U& z(f|O+!tKsF1HcxlMC*5MlQ32(o>LNF?Bw8b4g^67Q&d#+{a^q95o7=$AP_KYfH!Y| z(|#X7hB`pBVGuzh4`Lbs*bo7b2quC9Pjz%8t~TTEBz_PG#Kgn|fk629_#`GKQYe(_ z>S_jq!D6xae7+=v(AYRbkU<-s6uJt>K~ON31lG1z3IyN~Ac+QeN%q7eX-5pw&_8gR zl2WK76#!tTzMTQSrw9UX{{$7Xz(gk@kqDjx#KiC~xv|2@yY#w%`V<+F=zS4U9efH1 zk_6d6|5(vjhLeqJfFV%=fSGtvk#Qf7HFZ_c7JpJQbYc*-BdDm_%O8?pj^=&RaewT)Yh6eN%qN5LZ z>*@0~IMGI&YkmYZA1&rPs-#!RPnp~$o#GkP=Pk!V$Tb^-UQ01>JSM8C_~JeD z=NaNPuY`sD-k0`MP41@ATB6hDd*ZCv=!8jVn2WQXSF|7Z$73b zGz434a-*uROr}5x2it>Ft`!w(eeJ%g?pInGxqpcq=WNL=IhTJ-oAq}2O2P0$K}FdLKL95P%hu2RMgI(WTi`dY{OhWcBFhBy{fpl( z{73lvjSH_Wb&MF+>6G|*gMHdEDW*@61f3s7uiC~%4e!>BNgDBNchdds6_+w|T>N#| zGx|menQ_j=WwS_e!sd!=c~dxQ)gWRdF*6)>o7!(f0_^qy@U7}_KtqbCu_Z$TAZ`Db z{AU|t9-`wW-vX;SUSVGlbK|3%uX}*wj9C55i`{M~QCeGBdD~MM2TbM%m8)SnSgu?r za&|2xXua4XY_hnkr`>Eh^}Uu6bVj8;osG;QT3f6{U~al;Y0^p;osgsT!&S-&{P=zZ zL(KBbLVG`qUyXg37Vvs+;wtdgKrO?~b{!Gs=rz@J1ZlGKq9+n8<5s3W{>4$U3#4vPChBm#AiBis(cSb% z9WOc1%(&X7GB+N!_e-2JuWhqlwOvhFcHZuezPM-ExWsF0-;0GOJ!Tk{`n>0(jrmWo z@cfzC&y)$d%DN`ktYG?Pl(>3l9LJ|UF#l|Du1m6d6G^8(#eP-DKc(hWcYY2HKTt`v zvoo8tJ6VMqi>agvQa-iVP51oC%;%aG1#8W0CyWXjrM*#G$x!HM&H$_nR!`h4QI&ay zg-1q@$cR3-$|R^b!Pc<1HxSOlt`#(u4lhHKK>R3hChHB2G;O}rM2g+*$PXQ&t?4h$ zP%0~0JH7g?{OA`Q5J&3ljpZAl+fqxz2|LRtj()r|JwxgYlW zJ*J;9_1IR0Jo}{+OI9_1hOuU6Jhjbpih~I&pK)QI(;YFPI_$NDCyB{Yd+a3kj`QF? zIOgqM-7k7~5Srhp{m0C5T=1orZtBD!Eo2Y*?D%&w{Iiqc|Mt_|S0`xBAqeFg%&v_@ zDNmH6EH){J(<>i_(RhB>&l|`b2(LLz#vYp0Tt#GyQ*W7>1=ad;3t$efQ^l93Ba{}S z51i4^lSK$J2c}-6{Yt0C%C~Q1yrT@R`Tz7uz^V-!mA`pjcRs(SU-!Y{%MK2CE_^Eg z>i~K3nXLh?vSJ($hU9mxf7t5(@dTO~)6$pghGHjtwlA?fb=8SmADC5;=GoW_xiW99 z1|JNpSP9y+40_qGSK0G>qHb6F4)VxOL-pk=LF9Lh@HQ=cTllV_cj1dnmnh|zSVzjn z2kG(&#YT%rWM3${9@Qr}WG44BMOj*fR}m|V!r(gxQG>JC0NuWTqgsBmx-g{n`DdGyFDWY`QH{;`JVDo*jD#8b?w{-22VW{S96mU zl4xks##yYR4^B(1@qHqBK>zc}6cDm_YifsDd6Ua;7%x22cG{=Zjo$-dbtYZ6qC`fU z%O)@W(qPJr!hbFp*jo#UGhZs(MG`w}xvC;ELq@qWlU6)Qq-HX2_&D z?8d3e66Z@DZ(sva3pMX7E+d`u10S^_Lw1xd>1iT5$GF}4NM3-v{L7vjY0z)($YAN( zg6^_urTc4`nfc)g#ot=WR}8)s)b*g|7vjjn18kj<7hU>0qF46YtIal1V`d5uetU(V zTv>|te|a6nx%;h30wyFN9CgBVXr8suxiD}8?Vy=B5eiv9dH~(j{u-&%Y&gQFVyFY d>>qc+91!^BwiG6&Go|4>0Y8U4TVd;w^j|}Z)8YUC literal 0 HcmV?d00001 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 | -- 2.39.5