]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement hp_line_dof_identities in 2D and
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 10 Jul 2017 19:37:01 +0000 (21:37 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Wed, 12 Jul 2017 09:25:28 +0000 (11:25 +0200)
hp_quad_dof_identities in 3D for FE_FaceQ.

include/deal.II/fe/fe_face.h
source/fe/fe_face.cc
tests/hp/hp_line_dof_identities_faceq.cc [new file with mode: 0644]
tests/hp/hp_line_dof_identities_faceq.output [new file with mode: 0644]
tests/hp/hp_quad_dof_identities_faceq.cc [new file with mode: 0644]
tests/hp/hp_quad_dof_identities_faceq.output [new file with mode: 0644]
tests/hp/hp_vertex_dof_identities_faceq.cc [new file with mode: 0644]
tests/hp/hp_vertex_dof_identities_faceq.output [new file with mode: 0644]

index 3946e4a322acd4c47ce5539eb2b4dfe6f3bbbe44..4c1d42bc87267d4d3213e87182efda27ac7c7d80 100644 (file)
@@ -261,6 +261,50 @@ public:
    */
   virtual bool hp_constraints_are_implemented () const;
 
+  /**
+   * If, on a vertex, several finite elements are active, the hp code first
+   * assigns the degrees of freedom of each of these FEs different global
+   * indices. It then calls this function to find out which of them should get
+   * identical values, and consequently can receive the same global DoF index.
+   * This function therefore returns a list of identities between DoFs of the
+   * present finite element object with the DoFs of @p fe_other, which is a
+   * reference to a finite element object representing one of the other finite
+   * elements active on this particular vertex. The function computes which of
+   * the degrees of freedom of the two finite element objects are equivalent,
+   * both numbered between zero and the corresponding value of dofs_per_vertex
+   * of the two finite elements. The first index of each pair denotes one of
+   * the vertex dofs of the present element, whereas the second is the
+   * corresponding index of the other finite element.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_vertex_dof_identities (const FiniteElement<1, spacedim> &fe_other) const;
+
+  /**
+   * Same as hp_vertex_dof_indices(), except that the function treats degrees
+   * of freedom on lines.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_line_dof_identities (const FiniteElement<1, spacedim> &fe_other) const;
+
+  /**
+   * Same as hp_vertex_dof_indices(), except that the function treats degrees
+   * of freedom on quads.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_quad_dof_identities (const FiniteElement<1, spacedim> &fe_other) const;
+
   /**
    * Return whether this element dominates the one given as argument when they
    * meet at a common face, whether it is the other way around, whether
index eb1d5efbe666e77df2b4c16ec6896b88d82fa7a2..ed73ac5a99c7008e94b7b3347c7f1797ab591f1c 100644 (file)
@@ -268,10 +268,7 @@ std::vector<std::pair<unsigned int, unsigned int> >
 FE_FaceQ<dim, spacedim>::
 hp_vertex_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
 {
-  // this element is discontinuous, so by definition there can
-  // be no identities between its dofs and those of any neighbor
-  // (of whichever type the neighbor may be -- after all, we have
-  // no face dofs on this side to begin with)
+  // this element is always discontinuous at vertices
   return
     std::vector<std::pair<unsigned int, unsigned int> > ();
 }
@@ -281,14 +278,66 @@ hp_vertex_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) cons
 template <int dim, int spacedim>
 std::vector<std::pair<unsigned int, unsigned int> >
 FE_FaceQ<dim, spacedim>::
-hp_line_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+hp_line_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const
 {
-  // this element is discontinuous, so by definition there can
-  // be no identities between its dofs and those of any neighbor
-  // (of whichever type the neighbor may be -- after all, we have
-  // no face dofs on this side to begin with)
-  return
-    std::vector<std::pair<unsigned int, unsigned int> > ();
+  // this element is continuous only for the highest dimensional bounding object
+  if (dim !=2)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      // this is similar to the FE_Q_Base class
+      if (const FE_FaceQ<dim,spacedim> *fe_q_other
+          = dynamic_cast<const FE_FaceQ<dim,spacedim>*>(&fe_other))
+        {
+          // dofs are located along lines, so two dofs are identical if they are
+          // located at identical positions.
+          // Therefore, read the points in unit_support_points for the
+          // first coordinate direction. We take the lexicographic ordering of the
+          // points in the second direction (i.e., y-direction) since we know
+          // that the first p+1 dofs are located at the left (x=0) face.
+          const unsigned int p = this->degree;
+          const unsigned int q = fe_q_other->degree;
+
+          std::vector<std::pair<unsigned int, unsigned int> > identities;
+
+          const std::vector<unsigned int> &index_map_inverse=
+            this->poly_space.get_numbering_inverse();
+          const std::vector<unsigned int> &index_map_inverse_other=
+            fe_q_other->poly_space.get_numbering_inverse();
+
+          for (unsigned int i=0; i<p+1; ++i)
+            for (unsigned int j=0; j<q+1; ++j)
+              if (std::fabs(this->unit_support_points[index_map_inverse[i]][dim-1]-
+                            fe_q_other->unit_support_points[index_map_inverse_other[j]][dim-1])
+                  < 1e-14)
+                identities.emplace_back (i, j);
+
+          return identities;
+        }
+      else if (dynamic_cast<const FE_Nothing<dim>*>(&fe_other) != nullptr)
+        {
+          // the FE_Nothing has no degrees of freedom, so there are no
+          // equivalencies to be recorded
+          return std::vector<std::pair<unsigned int, unsigned int> > ();
+        }
+      else if (fe_other.dofs_per_face == 0)
+        {
+          // if the other element has no elements on faces at all,
+          // then it would be impossible to enforce any kind of
+          // continuity even if we knew exactly what kind of element
+          // we have -- simply because the other element declares
+          // that it is discontinuous because it has no DoFs on
+          // its faces. in that case, just state that we have no
+          // constraints to declare
+          return std::vector<std::pair<unsigned int, unsigned int> > ();
+        }
+      else
+        {
+          Assert (false, ExcNotImplemented());
+          return std::vector<std::pair<unsigned int, unsigned int> > ();
+        }
+    }
 }
 
 
@@ -296,14 +345,72 @@ hp_line_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
 template <int dim, int spacedim>
 std::vector<std::pair<unsigned int, unsigned int> >
 FE_FaceQ<dim, spacedim>::
-hp_quad_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+hp_quad_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const
 {
-  // this element is discontinuous, so by definition there can
-  // be no identities between its dofs and those of any neighbor
-  // (of whichever type the neighbor may be -- after all, we have
-  // no face dofs on this side to begin with)
-  return
-    std::vector<std::pair<unsigned int, unsigned int> > ();
+  // this element is continuous only for the highest dimensional bounding object
+  if (dim !=3)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      // this is similar to the FE_Q_Base class
+      if (const FE_FaceQ<dim,spacedim> *fe_q_other
+          = dynamic_cast<const FE_FaceQ<dim,spacedim>*>(&fe_other))
+        {
+          // this works exactly like the line case above, except that now we have
+          // to have two indices i1, i2 and j1, j2 to characterize the dofs on the
+          // face of each of the finite elements. since they are ordered
+          // lexicographically along the first line and we have a tensor product,
+          // the rest is rather straightforward
+          const unsigned int p = this->degree;
+          const unsigned int q = fe_q_other->degree;
+
+          std::vector<std::pair<unsigned int, unsigned int> > identities;
+
+          const std::vector<unsigned int> &index_map_inverse=
+            this->poly_space.get_numbering_inverse();
+          const std::vector<unsigned int> &index_map_inverse_other=
+            fe_q_other->poly_space.get_numbering_inverse();
+
+          std::vector<std::pair<unsigned int, unsigned int> > identities_1d;
+
+          for (unsigned int i=0; i<p+1; ++i)
+            for (unsigned int j=0; j<q+1; ++j)
+              if (std::fabs(this->unit_support_points[index_map_inverse[i]][dim-2]-
+                            fe_q_other->unit_support_points[index_map_inverse_other[j]][dim-2])
+                  < 1e-14)
+                identities_1d.emplace_back (i, j);
+
+          for (unsigned int n1=0; n1<identities_1d.size(); ++n1)
+            for (unsigned int n2=0; n2<identities_1d.size(); ++n2)
+              identities.emplace_back (identities_1d[n1].first*(p+1)+identities_1d[n2].first,
+                                       identities_1d[n1].second*(q+1)+identities_1d[n2].second);
+
+          return identities;
+        }
+      else if (dynamic_cast<const FE_Nothing<dim>*>(&fe_other) != nullptr)
+        {
+          // the FE_Nothing has no degrees of freedom, so there are no
+          // equivalencies to be recorded
+          return std::vector<std::pair<unsigned int, unsigned int> > ();
+        }
+      else if (fe_other.dofs_per_face == 0)
+        {
+          // if the other element has no elements on faces at all,
+          // then it would be impossible to enforce any kind of
+          // continuity even if we knew exactly what kind of element
+          // we have -- simply because the other element declares
+          // that it is discontinuous because it has no DoFs on
+          // its faces. in that case, just state that we have no
+          // constraints to declare
+          return std::vector<std::pair<unsigned int, unsigned int> > ();
+        }
+      else
+        {
+          Assert (false, ExcNotImplemented());
+          return std::vector<std::pair<unsigned int, unsigned int> > ();
+        }
+    }
 }
 
 
@@ -462,6 +569,40 @@ FE_FaceQ<1,spacedim>::hp_constraints_are_implemented () const
   return true;
 }
 
+template <int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceQ<1, spacedim>::
+hp_vertex_dof_identities (const FiniteElement<1, spacedim> &/*fe_other*/) const
+{
+  // this element is always discontinuous at vertices
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > (1,std::make_pair(0U,0U));
+}
+
+
+
+template <int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceQ<1, spacedim>::
+hp_line_dof_identities (const FiniteElement<1, spacedim> &fe_other) const
+{
+  // this element is continuous only for the highest dimensional bounding object
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceQ<1, spacedim>::
+hp_quad_dof_identities (const FiniteElement<1, spacedim> &fe_other) const
+{
+  // this element is continuous only for the highest dimensional bounding object
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
 
 
 template <int spacedim>
diff --git a/tests/hp/hp_line_dof_identities_faceq.cc b/tests/hp/hp_line_dof_identities_faceq.cc
new file mode 100644 (file)
index 0000000..1b52f23
--- /dev/null
@@ -0,0 +1,74 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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 FE_FaceQ::hp_line_dof_identities
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/fe/fe_face.h>
+#include <deal.II/base/quadrature_lib.h>
+
+#include <fstream>
+
+
+template <int dim>
+void test ()
+{
+  hp::FECollection<dim> fe_collection;
+  for (unsigned int i=1; i<8-dim; ++i)
+    fe_collection.push_back (FE_FaceQ<dim>(i));
+
+  for (unsigned int i=0; i<fe_collection.size(); ++i)
+    for (unsigned int j=0; j<fe_collection.size(); ++j)
+      {
+        const std::vector<std::pair<unsigned int, unsigned int> >
+        identities = fe_collection[i].hp_line_dof_identities (fe_collection[j]);
+
+        deallog << "Identities for "
+                << fe_collection[i].get_name() << " and "
+                << fe_collection[j].get_name() << ": "
+                << identities.size()
+                << std::endl;
+
+        for (unsigned int k=0; k<identities.size(); ++k)
+          {
+            Assert (identities[k].first < fe_collection[i].dofs_per_line,
+                    ExcInternalError());
+            Assert (identities[k].second < fe_collection[j].dofs_per_line,
+                    ExcInternalError());
+
+            deallog << identities[k].first << ' '
+                    << identities[k].second
+                    << std::endl;
+          }
+      }
+}
+
+
+
+int main ()
+{
+  initlog();
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/hp/hp_line_dof_identities_faceq.output b/tests/hp/hp_line_dof_identities_faceq.output
new file mode 100644 (file)
index 0000000..84959c8
--- /dev/null
@@ -0,0 +1,141 @@
+
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(1): 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(2): 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(3): 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(4): 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(5): 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(6): 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(1): 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(2): 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(3): 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(4): 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(5): 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(6): 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(1): 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(2): 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(3): 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(4): 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(5): 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(6): 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(1): 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(2): 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(3): 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(4): 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(5): 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(6): 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(1): 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(2): 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(3): 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(4): 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(5): 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(6): 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(1): 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(2): 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(3): 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(4): 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(5): 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(6): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(1): 2
+DEAL::0 0
+DEAL::1 1
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(2): 2
+DEAL::0 0
+DEAL::1 2
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(3): 2
+DEAL::0 0
+DEAL::1 3
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(4): 2
+DEAL::0 0
+DEAL::1 4
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(5): 2
+DEAL::0 0
+DEAL::1 5
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(1): 2
+DEAL::0 0
+DEAL::2 1
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(2): 3
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(3): 2
+DEAL::0 0
+DEAL::2 3
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(4): 3
+DEAL::0 0
+DEAL::1 2
+DEAL::2 4
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(5): 2
+DEAL::0 0
+DEAL::2 5
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(1): 2
+DEAL::0 0
+DEAL::3 1
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(2): 2
+DEAL::0 0
+DEAL::3 2
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(3): 4
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(4): 2
+DEAL::0 0
+DEAL::3 4
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(5): 2
+DEAL::0 0
+DEAL::3 5
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(1): 2
+DEAL::0 0
+DEAL::4 1
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(2): 3
+DEAL::0 0
+DEAL::2 1
+DEAL::4 2
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(3): 2
+DEAL::0 0
+DEAL::4 3
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(4): 5
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(5): 2
+DEAL::0 0
+DEAL::4 5
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(1): 2
+DEAL::0 0
+DEAL::5 1
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(2): 2
+DEAL::0 0
+DEAL::5 2
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(3): 2
+DEAL::0 0
+DEAL::5 3
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(4): 2
+DEAL::0 0
+DEAL::5 4
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(5): 6
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 5
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(4): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(4): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(4): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(4): 0
+DEAL::OK
diff --git a/tests/hp/hp_quad_dof_identities_faceq.cc b/tests/hp/hp_quad_dof_identities_faceq.cc
new file mode 100644 (file)
index 0000000..1cc6920
--- /dev/null
@@ -0,0 +1,72 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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 FE_FaceQ::hp_quad_dof_identities
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/fe/fe_face.h>
+
+#include <fstream>
+
+
+template <int dim>
+void test ()
+{
+  hp::FECollection<dim> fe_collection;
+  for (unsigned int i=1; i<8-dim; ++i)
+    fe_collection.push_back (FE_FaceQ<dim>(i));
+
+  for (unsigned int i=0; i<fe_collection.size(); ++i)
+    for (unsigned int j=0; j<fe_collection.size(); ++j)
+      {
+        const std::vector<std::pair<unsigned int, unsigned int> >
+        identities = fe_collection[i].hp_quad_dof_identities (fe_collection[j]);
+
+        deallog << "Identities for "
+                << fe_collection[i].get_name() << " and "
+                << fe_collection[j].get_name() << ": "
+                << identities.size()
+                << std::endl;
+
+        for (unsigned int k=0; k<identities.size(); ++k)
+          {
+            Assert (identities[k].first < fe_collection[i].dofs_per_quad,
+                    ExcInternalError());
+            Assert (identities[k].second < fe_collection[j].dofs_per_quad,
+                    ExcInternalError());
+
+            deallog << identities[k].first << ' '
+                    << identities[k].second
+                    << std::endl;
+          }
+      }
+}
+
+
+
+int main ()
+{
+  initlog();
+
+  test<2> ();
+  test<3> ();
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/hp/hp_quad_dof_identities_faceq.output b/tests/hp/hp_quad_dof_identities_faceq.output
new file mode 100644 (file)
index 0000000..ceeabbe
--- /dev/null
@@ -0,0 +1,155 @@
+
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(1): 4
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(2): 4
+DEAL::0 0
+DEAL::1 2
+DEAL::2 6
+DEAL::3 8
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(3): 4
+DEAL::0 0
+DEAL::1 3
+DEAL::2 12
+DEAL::3 15
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(4): 4
+DEAL::0 0
+DEAL::1 4
+DEAL::2 20
+DEAL::3 24
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(1): 4
+DEAL::0 0
+DEAL::2 1
+DEAL::6 2
+DEAL::8 3
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(2): 9
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 5
+DEAL::6 6
+DEAL::7 7
+DEAL::8 8
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(3): 4
+DEAL::0 0
+DEAL::2 3
+DEAL::6 12
+DEAL::8 15
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(4): 9
+DEAL::0 0
+DEAL::1 2
+DEAL::2 4
+DEAL::3 10
+DEAL::4 12
+DEAL::5 14
+DEAL::6 20
+DEAL::7 22
+DEAL::8 24
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(1): 4
+DEAL::0 0
+DEAL::3 1
+DEAL::12 2
+DEAL::15 3
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(2): 4
+DEAL::0 0
+DEAL::3 2
+DEAL::12 6
+DEAL::15 8
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(3): 16
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 5
+DEAL::6 6
+DEAL::7 7
+DEAL::8 8
+DEAL::9 9
+DEAL::10 10
+DEAL::11 11
+DEAL::12 12
+DEAL::13 13
+DEAL::14 14
+DEAL::15 15
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(4): 4
+DEAL::0 0
+DEAL::3 4
+DEAL::12 20
+DEAL::15 24
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(1): 4
+DEAL::0 0
+DEAL::4 1
+DEAL::20 2
+DEAL::24 3
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(2): 9
+DEAL::0 0
+DEAL::2 1
+DEAL::4 2
+DEAL::10 3
+DEAL::12 4
+DEAL::14 5
+DEAL::20 6
+DEAL::22 7
+DEAL::24 8
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(3): 4
+DEAL::0 0
+DEAL::4 3
+DEAL::20 12
+DEAL::24 15
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(4): 25
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 5
+DEAL::6 6
+DEAL::7 7
+DEAL::8 8
+DEAL::9 9
+DEAL::10 10
+DEAL::11 11
+DEAL::12 12
+DEAL::13 13
+DEAL::14 14
+DEAL::15 15
+DEAL::16 16
+DEAL::17 17
+DEAL::18 18
+DEAL::19 19
+DEAL::20 20
+DEAL::21 21
+DEAL::22 22
+DEAL::23 23
+DEAL::24 24
+DEAL::OK
diff --git a/tests/hp/hp_vertex_dof_identities_faceq.cc b/tests/hp/hp_vertex_dof_identities_faceq.cc
new file mode 100644 (file)
index 0000000..1051e27
--- /dev/null
@@ -0,0 +1,73 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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 FE_FaceQ::hp_vertex_dof_identities
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/fe/fe_face.h>
+
+#include <fstream>
+
+
+template <int dim>
+void test ()
+{
+  hp::FECollection<dim> fe_collection;
+  for (unsigned int i=1; i<8-dim; ++i)
+    fe_collection.push_back (FE_FaceQ<dim>(i));
+
+  for (unsigned int i=0; i<fe_collection.size(); ++i)
+    for (unsigned int j=0; j<fe_collection.size(); ++j)
+      {
+        const std::vector<std::pair<unsigned int, unsigned int> >
+        identities = fe_collection[i].hp_vertex_dof_identities (fe_collection[j]);
+
+        deallog << "Identities for "
+                << fe_collection[i].get_name() << " and "
+                << fe_collection[j].get_name() << ": "
+                << identities.size()
+                << std::endl;
+
+        for (unsigned int k=0; k<identities.size(); ++k)
+          {
+            Assert (identities[k].first < fe_collection[i].dofs_per_vertex,
+                    ExcInternalError());
+            Assert (identities[k].second < fe_collection[j].dofs_per_vertex,
+                    ExcInternalError());
+
+            deallog << identities[k].first << ' '
+                    << identities[k].second
+                    << std::endl;
+          }
+      }
+}
+
+
+
+int main ()
+{
+  initlog();
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/hp/hp_vertex_dof_identities_faceq.output b/tests/hp/hp_vertex_dof_identities_faceq.output
new file mode 100644 (file)
index 0000000..d4a29e0
--- /dev/null
@@ -0,0 +1,115 @@
+
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(1): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(2): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(3): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(4): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(5): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(1) and FE_FaceQ<1>(6): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(1): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(2): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(3): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(4): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(5): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(2) and FE_FaceQ<1>(6): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(1): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(2): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(3): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(4): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(5): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(3) and FE_FaceQ<1>(6): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(1): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(2): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(3): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(4): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(5): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(4) and FE_FaceQ<1>(6): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(1): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(2): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(3): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(4): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(5): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(5) and FE_FaceQ<1>(6): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(1): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(2): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(3): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(4): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(5): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<1>(6) and FE_FaceQ<1>(6): 1
+DEAL::0 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(1) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(2) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(3) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(4) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(1): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(2): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(3): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(4): 0
+DEAL::Identities for FE_FaceQ<2>(5) and FE_FaceQ<2>(5): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(1) and FE_FaceQ<3>(4): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(2) and FE_FaceQ<3>(4): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(3) and FE_FaceQ<3>(4): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(1): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(2): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(3): 0
+DEAL::Identities for FE_FaceQ<3>(4) and FE_FaceQ<3>(4): 0
+DEAL::OK

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.