]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Allow querying whether a DoF is located on a vertex, line, etc. 153/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 17 Sep 2014 14:42:44 +0000 (09:42 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 30 Sep 2014 23:13:53 +0000 (18:13 -0500)
doc/news/changes.h
include/deal.II/base/geometry_info.h
include/deal.II/fe/fe.h
tests/fe/get_dof_association.cc [new file with mode: 0644]
tests/fe/get_dof_association.output [new file with mode: 0644]

index 4b6a158742697e70f55383d2c876124a5f63a772..0665d082ae3f867faf50b7fa93c455371c6074fb 100644 (file)
@@ -307,6 +307,13 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> New: The new function FiniteElement::get_associated_geometry_primitive() allows to
+  query whether a given degree of freedom is associated with a vertex, line,
+  quad, or hex.
+  <br>
+  (Wolfgang Bangerth, 2014/09/26)
+  </li>
+
   <li> Fixed: The vector and array versions of Utilities::MPI::sum() and
   Utilities::MPI::max() produced segmentation faults with some MPI implementations
   if the input and output arguments were the same. This is now fixed.
index c1eb5547f3a4ba1f5da38374c54600a391d7808d..2694688cd6014379b3dfe985208aece74a02cfca 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 1998 - 2013 by the deal.II authors
+// Copyright (C) 1998 - 2014 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
 DEAL_II_NAMESPACE_OPEN
 
 
+/**
+ * A class that can represent the kinds of objects a triangulation is
+ * made up of: vertices, lines, quads and hexes.
+ *
+ * The class is rather primitive: it only stores a single integer
+ * that represents the dimensionality of the object represented.
+ * In other words, this class is useful primarily as a way to
+ * pass around an object whose data type explains what it does
+ * (unlike just passing around an integer), and for providing
+ * symbolic names for these objects such as GeometryPrimitive::vertex
+ * instead of an integer zero.
+ *
+ * Since the ability to identify such objects with the integral
+ * dimension of the object represented, this class provides
+ * conversion operators to and from unsigned integers.
+ *
+ * @author Wolfgang Bangerth, 2014
+ */
+class GeometryPrimitive
+{
+public:
+  /**
+   * An enumeration providing symbolic names for the
+   * objects that can be represented by this class.
+   * The numeric values of these symbolic names equal
+   * the geometric dimensionality of the represented
+   * objects to make conversion from and to integer
+   * variables simpler.
+   */
+  enum Object
+  {
+    vertex = 0,
+    line   = 1,
+    quad   = 2,
+    hex    = 3
+  };
+
+  /**
+   * Constructor. Initialize the object with the
+   * given argument representing a vertex, line, etc.
+   */
+  GeometryPrimitive (const Object object);
+
+  /**
+   * Constructor. Initialize the object with an
+   * integer that should represent the dimensionality of
+   * the geometric object in question. This will usually be
+   * a number between zero (a vertex) and three (a hexahedron).
+   */
+  GeometryPrimitive (const unsigned int object_dimension);
+
+  /**
+   * Return the integral dimension of the object
+   * currently represented, i.e. zero for a vertex,
+   * one for a line, etc.
+   */
+  operator unsigned int () const;
+
+private:
+  /**
+   * The object currently represented.
+   */
+  Object object;
+};
+
+
+
 
 /**
  * A class that provides possible choices for isotropic and
@@ -2316,6 +2383,30 @@ d_linear_shape_function_gradient (const Point<3> &xi,
 
 /* -------------- inline functions ------------- */
 
+
+inline
+GeometryPrimitive::GeometryPrimitive (const Object object)
+  :
+  object (object)
+{}
+
+
+
+inline
+GeometryPrimitive::GeometryPrimitive (const unsigned int object_dimension)
+  :
+  object (static_cast<Object>(object_dimension))
+{}
+
+
+inline
+GeometryPrimitive::operator unsigned int () const
+{
+  return static_cast<unsigned int>(object);
+}
+
+
+
 namespace internal
 {
 
index c664461caf37554a041a207fcdb60546dc88dfcd..a6e658469f7f3aef9848763bb4bec8969a793f97 100644 (file)
@@ -1499,6 +1499,56 @@ public:
   bool
   has_generalized_face_support_points () const;
 
+  /**
+   * For a given degree of freedom, return whether it is logically
+   * associated with a vertex, line, quad or hex.
+   *
+   * For instance, for continuous finite elements this coincides with
+   * the lowest dimensional object the support point of the degree of
+   * freedom lies on. To give an example, for $Q_1$ elements in 3d,
+   * every degree of freedom is defined by a shape function that we
+   * get by interpolating using support points that lie on the
+   * vertices of the cell. The support of these points of course
+   * extends to all edges connected to this vertex, as well
+   * as the adjacent faces and the cell interior, but we say that
+   * logically the degree of freedom is associated with the vertex as
+   * this is the lowest-dimensional object it is associated with.
+   * Likewise, for $Q_2$ elements in 3d, the degrees of freedom
+   * with support points at edge midpoints would yield a value of
+   * GeometryPrimitive::line from this function, whereas those on
+   * the centers of faces in 3d would return GeometryPrimitive::quad.
+   *
+   * To make this more formal, the kind of object returned by
+   * this function represents the object so that the support of the
+   * shape function corresponding to the degree of freedom, (i.e., that
+   * part of the domain where the function "lives") is the union
+   * of all of the cells sharing this object. To return to the example
+   * above, for $Q_2$ in 3d, the shape function with support point at
+   * an edge midpoint has support on all cells that share the edge and
+   * not only the cells that share the adjacent faces, and consequently
+   * the function will return GeometryPrimitive::line.
+   *
+   * On the other hand, for discontinuous elements of type $DGQ_2$,
+   * a degree of freedom associated with an interpolation polynomial
+   * that has its support point physically located at a line bounding
+   * a cell, but is nonzero only on one cell. Consequently, it is
+   * logically associated with the interior of that cell (i.e., with a
+   * GeometryPrimitive::quad in 2d and a GeometryPrimitive::hex in 3d).
+   *
+   * @param[in] cell_dof_index The index of a shape function or
+   *   degree of freedom. This index must be in the range
+   *   <code>[0,dofs_per_cell)</code>.
+   *
+   * @note The integer value of the object returned by this function
+   *   equals the dimensionality of the object it describes, and can
+   *   consequently be used in generic programming paradigms. For
+   *   example, if a degree of freedom is associated with a vertex,
+   *   then this function returns GeometryPrimitive::vertex, which has
+   *   a numeric value of zero (the dimensionality of a vertex).
+   */
+  GeometryPrimitive
+  get_associated_geometry_primitive (const unsigned int cell_dof_index) const;
+
   /**
    * Interpolate a set of scalar values, computed in the generalized support
    * points.
@@ -2283,6 +2333,28 @@ FiniteElement<dim,spacedim>::is_primitive (const unsigned int i) const
 
 
 
+template <int dim, int spacedim>
+inline
+GeometryPrimitive
+FiniteElement<dim,spacedim>::get_associated_geometry_primitive (const unsigned int cell_dof_index) const
+{
+  Assert (cell_dof_index < this->dofs_per_cell,
+          ExcIndexRange (cell_dof_index, 0, this->dofs_per_cell));
+
+  // just go through the usual cases, taking into account how DoFs
+  // are enumerated on the reference cell
+  if (cell_dof_index < this->first_line_index)
+    return GeometryPrimitive::vertex;
+  else if (cell_dof_index < this->first_quad_index)
+    return GeometryPrimitive::line;
+  else if (cell_dof_index < this->first_hex_index)
+    return GeometryPrimitive::quad;
+  else
+    return GeometryPrimitive::hex;
+}
+
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
diff --git a/tests/fe/get_dof_association.cc b/tests/fe/get_dof_association.cc
new file mode 100644 (file)
index 0000000..df67185
--- /dev/null
@@ -0,0 +1,72 @@
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 1998 - 2014 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_dof_association() with a couple of elements
+
+
+#include "../tests.h"
+#include <iostream>
+#include <fstream>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_nedelec.h>
+
+
+template <int dim>
+void test (const FiniteElement<dim> &fe)
+{
+  deallog << fe.get_name() << std::endl;
+
+  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
+    {
+      switch (fe.get_associated_geometry_primitive(i))
+        {
+        case GeometryPrimitive::vertex:
+          deallog << 'v';
+          break;
+        case GeometryPrimitive::line:
+          deallog << 'l';
+          break;
+        case GeometryPrimitive::quad:
+          deallog << 'q';
+          break;
+        case GeometryPrimitive::hex:
+          deallog << 'h';
+          break;
+        default:
+          Assert (false, ExcInternalError());
+        }
+    }
+  deallog << std::endl;
+}
+
+int main()
+{
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test(FE_Q<1>(2));
+  test(FE_Q<2>(2));
+  test(FE_Q<3>(2));
+
+  test(FE_Nedelec<2>(0));
+  test(FE_Nedelec<3>(0));
+
+  test(FE_Nedelec<2>(1));
+  test(FE_Nedelec<3>(1));
+}
diff --git a/tests/fe/get_dof_association.output b/tests/fe/get_dof_association.output
new file mode 100644 (file)
index 0000000..01f8372
--- /dev/null
@@ -0,0 +1,15 @@
+
+DEAL::FE_Q<1>(2)
+DEAL::vvl
+DEAL::FE_Q<2>(2)
+DEAL::vvvvllllq
+DEAL::FE_Q<3>(2)
+DEAL::vvvvvvvvllllllllllllqqqqqqh
+DEAL::FE_Nedelec<2>(0)
+DEAL::llll
+DEAL::FE_Nedelec<3>(0)
+DEAL::llllllllllll
+DEAL::FE_Nedelec<2>(1)
+DEAL::llllllllqqqq
+DEAL::FE_Nedelec<3>(1)
+DEAL::llllllllllllllllllllllllqqqqqqqqqqqqqqqqqqqqqqqqhhhhhh

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.