]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clean up FE_Nothing::get_name(). 11437/head
authorDavid Wells <drwells@email.unc.edu>
Thu, 24 Dec 2020 02:14:45 +0000 (21:14 -0500)
committerDavid Wells <drwells@email.unc.edu>
Wed, 30 Dec 2020 23:38:36 +0000 (18:38 -0500)
include/deal.II/fe/fe_nothing.h
include/deal.II/grid/reference_cell.h
source/fe/fe_nothing.cc
tests/fe/get_fe_by_name_01.output
tests/hp/fe_nothing_22.cc
tests/hp/fe_nothing_22.output

index b3b83ec8d37735df79b4db94e88af8c45e6220a1..dd779d1ee311bbab8e11e38f4f07234d0dd994d3 100644 (file)
@@ -158,8 +158,17 @@ public:
   clone() const override;
 
   /**
-   * Return a string that uniquely identifies a finite element. In this case
-   * it is <code>FE_Nothing@<dim@></code>.
+   * Return a string that uniquely identifies a finite element. The name is
+   * <tt>FE_Nothing@<dim,spacedim@>(type, n_components, dominating)</tt> where
+   * <tt>dim</tt>, <tt>spacedim</tt>, <tt>type</tt>, and <tt>n_components</tt>
+   * are all specified by the constructor or type signature with the following
+   * exceptions:
+   * <ol>
+   *   <li>If <tt>spacedim == dim</tt> then that field is not printed.</li>
+   *   <li>If <tt>type</tt> is a hypercube then that field is not printed.</li>
+   *   <li>If <tt>n_components == 1</tt> then that field is not printed.</li>
+   *   <li>If <tt>dominate == false</tt> then that field is not printed.</li>
+   * </ol>
    */
   virtual std::string
   get_name() const override;
@@ -323,16 +332,6 @@ public:
   bool
   is_dominating() const;
 
-  /**
-   * Comparison operator. In addition to the fields already checked by
-   * FiniteElement::operator==(), this operator also checks for equality
-   * of the arguments passed to the constructors of the current object
-   * as well as the object against which the comparison is done (which
-   * for this purpose obviously also needs to be of type FE_Nothing).
-   */
-  virtual bool
-  operator==(const FiniteElement<dim, spacedim> &fe) const override;
-
 private:
   /**
    * If true, this element will dominate any other apart from itself in
index 2bfa0943b121fe729a3fccaa587c60a57efb49c8..2cd1399533574fec61e6c037d48ceab227188353 100644 (file)
@@ -70,6 +70,39 @@ namespace ReferenceCell
       }
   }
 
+  /**
+   * Convert the given reference cell type to a string.
+   */
+  inline std::string
+  to_string(const Type &type)
+  {
+    switch (type)
+      {
+        case Type::Vertex:
+          return "Vertex";
+        case Type::Line:
+          return "Line";
+        case Type::Tri:
+          return "Tri";
+        case Type::Quad:
+          return "Quad";
+        case Type::Tet:
+          return "Tet";
+        case Type::Pyramid:
+          return "Pyramid";
+        case Type::Wedge:
+          return "Wedge";
+        case Type::Hex:
+          return "Hex";
+        case Type::Invalid:
+          return "Invalid";
+        default:
+          Assert(false, ExcNotImplemented());
+      }
+
+    return "Invalid";
+  }
+
   /**
    * Return the correct simplex reference cell type for the given dimension
    * @p dim.
index 5cba50338de7598c1db459078479e92aeec95a88..44009a9e489c0979d77b6d34c8387d3620bc737a 100644 (file)
@@ -72,16 +72,25 @@ std::string
 FE_Nothing<dim, spacedim>::get_name() const
 {
   std::ostringstream namebuf;
-  namebuf << "FE_Nothing<" << dim << ">(";
+  namebuf << "FE_Nothing<" << Utilities::dim_string(dim, spacedim) << ">(";
+
+  std::vector<std::string> name_components;
+  if (this->reference_cell_type() != ReferenceCell::get_hypercube(dim))
+    name_components.push_back(
+      ReferenceCell::to_string(this->reference_cell_type()));
   if (this->n_components() > 1)
+    name_components.push_back(std::to_string(this->n_components()));
+  if (dominate)
+    name_components.emplace_back("dominating");
+
+  for (const std::string &comp : name_components)
     {
-      namebuf << this->n_components();
-      if (dominate)
-        namebuf << ", dominating";
+      namebuf << comp;
+      if (comp != name_components.back())
+        namebuf << ", ";
     }
-  else if (dominate)
-    namebuf << "dominating";
   namebuf << ")";
+
   return namebuf.str();
 }
 
@@ -199,29 +208,6 @@ FE_Nothing<dim, spacedim>::is_dominating() const
 
 
 
-template <int dim, int spacedim>
-bool
-FE_Nothing<dim, spacedim>::
-operator==(const FiniteElement<dim, spacedim> &f) const
-{
-  // Compare fields stored in the base class
-  if (!(this->FiniteElement<dim, spacedim>::operator==(f)))
-    return false;
-
-  // Then make sure the other object is really of type FE_Nothing,
-  // and compare the data that has been passed to both objects'
-  // constructors.
-  if (const FE_Nothing<dim, spacedim> *fe_nothing =
-        dynamic_cast<const FE_Nothing<dim, spacedim> *>(&f))
-    return ((dominate == fe_nothing->dominate) &&
-            (this->components == fe_nothing->components) &&
-            (this->reference_cell_type() == fe_nothing->reference_cell_type()));
-  else
-    return false;
-}
-
-
-
 template <int dim, int spacedim>
 FiniteElementDomination::Domination
 FE_Nothing<dim, spacedim>::compare_for_domination(
index 096155d2c03bdb8927e6d98c20343024c80ad944..f2249c8265d4b4a9c26390b951ce69b729ee78e2 100644 (file)
@@ -271,16 +271,16 @@ DEAL::Generated :
 DEAL::FE_Nothing<1>()
 DEAL::Read FE_Nothing()
 DEAL::Generated :
-DEAL::FE_Nothing<1>()
+DEAL::FE_Nothing<1,2>()
 DEAL::Read FE_Nothing()
 DEAL::Generated :
 DEAL::FE_Nothing<2>()
 DEAL::Read FE_Nothing()
 DEAL::Generated :
-DEAL::FE_Nothing<1>()
+DEAL::FE_Nothing<1,3>()
 DEAL::Read FE_Nothing()
 DEAL::Generated :
-DEAL::FE_Nothing<2>()
+DEAL::FE_Nothing<2,3>()
 DEAL::Read FE_Nothing()
 DEAL::Generated :
 DEAL::FE_Nothing<3>()
index 7bc1cc828d9ff27dba363d2002dc60925a23cd1a..eb24c88f85af503665363076c5dcca53cf73f434 100644 (file)
 
 
 
-// Test FE_Nothing::operator==()
-
-
-#include <deal.II/base/function.h>
-#include <deal.II/base/quadrature_lib.h>
-
-#include <deal.II/dofs/dof_accessor.h>
-#include <deal.II/dofs/dof_handler.h>
-#include <deal.II/dofs/dof_tools.h>
+// Test FE_Nothing::operator==(). The base clase operator should suffice
 
 #include <deal.II/fe/fe_nothing.h>
-#include <deal.II/fe/fe_q.h>
-#include <deal.II/fe/fe_system.h>
-
-#include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/grid_refinement.h>
-#include <deal.II/grid/tria.h>
-#include <deal.II/grid/tria_accessor.h>
-#include <deal.II/grid/tria_iterator.h>
 
-#include <deal.II/numerics/vector_tools.h>
+#include <deal.II/grid/reference_cell.h>
 
 #include "../tests.h"
 
@@ -45,6 +29,7 @@ template <int dim>
 void
 test()
 {
+  deallog << "dim = " << dim << std::endl;
   deallog << std::boolalpha;
   deallog << (FE_Nothing<dim>(1) == FE_Nothing<dim>(1, false)) << std::endl;
   deallog << (FE_Nothing<dim>(1) == FE_Nothing<dim>(2)) << std::endl;
@@ -52,6 +37,30 @@ test()
           << std::endl;
   deallog << (FE_Nothing<dim>(1, true) == FE_Nothing<dim>(2, true))
           << std::endl;
+  if (dim == 2)
+    {
+      deallog << (FE_Nothing<dim>(ReferenceCell::Type::Quad, 2, true) ==
+                  FE_Nothing<dim>(2, true))
+              << std::endl;
+      deallog << (FE_Nothing<dim>(ReferenceCell::Type::Tri, 2, true) ==
+                  FE_Nothing<dim>(2, true))
+              << std::endl;
+    }
+  if (dim == 3)
+    {
+      deallog << (FE_Nothing<dim>(ReferenceCell::Type::Hex, 2, true) ==
+                  FE_Nothing<dim>(2, true))
+              << std::endl;
+      deallog << (FE_Nothing<dim>(ReferenceCell::Type::Tet, 2, true) ==
+                  FE_Nothing<dim>(2, true))
+              << std::endl;
+      deallog << (FE_Nothing<dim>(ReferenceCell::Type::Wedge, 1, false) ==
+                  FE_Nothing<dim>(ReferenceCell::Type::Pyramid, 1, false))
+              << std::endl;
+      deallog << (FE_Nothing<dim>(ReferenceCell::Type::Wedge, 3) ==
+                  FE_Nothing<dim>(3))
+              << std::endl;
+    }
 }
 
 
index 9ce031bab882a2b122e182b201648ddb556ae354..17d12fbb372acfe94df7161a603fcb6a35332fbf 100644 (file)
@@ -1,8 +1,17 @@
 
+DEAL::dim = 1
 DEAL::true
 DEAL::false
 DEAL::false
 DEAL::false
+DEAL::dim = 2
+DEAL::true
+DEAL::false
+DEAL::false
+DEAL::false
+DEAL::true
+DEAL::false
+DEAL::dim = 3
 DEAL::true
 DEAL::false
 DEAL::false

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.