Added "How do I get the degree of freedom indices at vertices"
authorJean-Paul Pelteret <jppelteret@gmail.com>
Mon, 25 Feb 2019 08:24:51 +0000 (09:24 +0100)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Mon, 25 Feb 2019 08:24:51 +0000 (09:24 +0100)
Frequently-Asked-Questions.md

index bc74b0b..458f4f5 100644 (file)
@@ -50,6 +50,7 @@ This page collects a few answers to questions that have frequently been asked ab
     * [Questions about specific behavior of parts of deal.II](#questions-about-specific-behavior-of-parts-of-dealii)
       * [How do I create the mesh for my problem?](#how-do-i-create-the-mesh-for-my-problem)
       * [How do I describe complex boundaries?](#how-do-i-describe-complex-boundaries)
+      * [How do I get the degree of freedom indices at vertices?](#how-do-i-get-the-degree-of-freedom-indices-at-vertices) 
       * [I am using discontinuous Lagrange elements (FE_DGQ) but they don't seem to have vertex degrees of freedom!?](#i-am-using-discontinuous-lagrange-elements-fe_dgq-but-they-dont-seem-to-have-vertex-degrees-of-freedom)
       * [How do I access values of discontinuous elements at vertices?](#how-do-i-access-values-of-discontinuous-elements-at-vertices)
       * [Does deal.II support anisotropic finite element shape functions?](#does-dealii-support-anisotropic-finite-element-shape-functions)
@@ -1433,6 +1434,35 @@ the interior lie. The step-53 tutorial program explains how this is done
 for a realistic example.
 
 
+### How do I get the degree of freedom indices at vertices?
+
+For simple cases (for example, where only `FE_Q` elements are used) you could 
+use the `cell->vertex_dof_index()` function. This would mean that you'd require 
+a single grid traversal to extract each DoF value  associated with support 
+points corresponding to the vertices. The code to do this would look something 
+like this: 
+```cpp
+auto cell = dof_handler.begin_active ();
+const auto endc = dof_handler.end ();
+for (; cell != endc; ++cell)
+{
+  for (unsigned int vertex = 0; vertex < GeometryInfo<dim>::vertices_per_cell; 
+       ++vertex)
+  {
+    for (unsigned int component=0; component<fe.n_components(); ++component)
+      // Index associated with the given component at the cell local vertex
+      const unsigned int idx = cell->vertex_dof_index(vertex,component);
+  }
+}  
+```
+If this is something that you might repeat often then you could use the `GridTools::find_cells_adjacent_to_vertex()` function to cache this association 
+once up front for a single vertex, or the `GridTools::vertex_to_cell_map()` function 
+to do the same for all vertices. 
+
+For the most general case (e.g., when using non-primitive finite elements), you 
+might need to use `DoFTools::map_dofs_to_support_points()` function and then find
+which support points at located at the vertex position that you're interested in.
+
 ### I am using discontinuous Lagrange elements (`FE_DGQ`) but they don't seem to have vertex degrees of freedom!?
 
 Indeed. And here's the reason: a vertex is an entity that is shared between

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.