]> https://gitweb.dealii.org/ - dealii.git/commitdiff
hp for RTNodal and FENothing 4996/head
authorEldar Khattatov <eld.khattatov@gmail.com>
Thu, 31 Aug 2017 21:46:15 +0000 (17:46 -0400)
committerEldar Khattatov <eld.khattatov@gmail.com>
Thu, 31 Aug 2017 21:46:15 +0000 (17:46 -0400)
source/fe/fe_raviart_thomas_nodal.cc
tests/hp/face_domination_04.cc [new file with mode: 0644]
tests/hp/face_domination_04.output [new file with mode: 0644]

index ff145a40b8478a764fc65aa84cf8868d19a84d61..b6583308aec94aada173b73053137d450bdb261c 100644 (file)
@@ -22,6 +22,7 @@
 #include <deal.II/dofs/dof_accessor.h>
 #include <deal.II/fe/fe.h>
 #include <deal.II/fe/mapping.h>
+#include <deal.II/fe/fe_nothing.h>
 #include <deal.II/fe/fe_raviart_thomas.h>
 #include <deal.II/fe/fe_values.h>
 #include <deal.II/fe/fe_tools.h>
@@ -347,11 +348,13 @@ FE_RaviartThomasNodal<dim>::hp_vertex_dof_identities (
 {
   // we can presently only compute these
   // identities if both FEs are
-  // FE_RaviartThomasNodals. in that case, no
-  // dofs are assigned on the vertex, so we
-  // shouldn't be getting here at all.
+  // FE_RaviartThomasNodals or the other is FE_Nothing.
+  // In either case, no dofs are assigned on the vertex,
+  // so we shouldn't be getting here at all.
   if (dynamic_cast<const FE_RaviartThomasNodal<dim>*>(&fe_other)!=nullptr)
     return std::vector<std::pair<unsigned int, unsigned int> > ();
+  else if (dynamic_cast<const FE_Nothing<dim>*>(&fe_other) != nullptr)
+    return std::vector<std::pair<unsigned int, unsigned int> > ();
   else
     {
       Assert(false, ExcNotImplemented());
@@ -368,7 +371,8 @@ hp_line_dof_identities (const FiniteElement<dim> &fe_other) const
 {
   // we can presently only compute
   // these identities if both FEs are
-  // FE_RaviartThomasNodals
+  // FE_RaviartThomasNodals or if the other
+  // one is FE_Nothing
   if (const FE_RaviartThomasNodal<dim> *fe_q_other
       = dynamic_cast<const FE_RaviartThomasNodal<dim>*>(&fe_other))
     {
@@ -410,6 +414,12 @@ hp_line_dof_identities (const FiniteElement<dim> &fe_other) const
 
       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
     {
       Assert (false, ExcNotImplemented());
@@ -425,7 +435,8 @@ FE_RaviartThomasNodal<dim>::hp_quad_dof_identities (
 {
   // we can presently only compute
   // these identities if both FEs are
-  // FE_RaviartThomasNodals
+  // FE_RaviartThomasNodals or if the other
+  // one is FE_Nothing
   if (const FE_RaviartThomasNodal<dim> *fe_q_other
       = dynamic_cast<const FE_RaviartThomasNodal<dim>*>(&fe_other))
     {
@@ -450,6 +461,12 @@ FE_RaviartThomasNodal<dim>::hp_quad_dof_identities (
 
       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
     {
       Assert (false, ExcNotImplemented());
@@ -473,6 +490,21 @@ FE_RaviartThomasNodal<dim>::compare_for_face_domination (
       else
         return FiniteElementDomination::other_element_dominates;
     }
+  else if (const FE_Nothing<dim> *fe_q_other
+           = dynamic_cast<const FE_Nothing<dim>*>(&fe_other))
+    {
+      if (fe_q_other->is_dominating())
+        {
+          return FiniteElementDomination::other_element_dominates;
+        }
+      else
+        {
+          // FE_Nothing has no degrees of freedom and is typically
+          // used in a context where there are no continuity
+          // requirements along the interface
+          return FiniteElementDomination::no_requirements;
+        }
+    }
 
   Assert (false, ExcNotImplemented());
   return FiniteElementDomination::neither_element_dominates;
diff --git a/tests/hp/face_domination_04.cc b/tests/hp/face_domination_04.cc
new file mode 100644 (file)
index 0000000..f49f1d7
--- /dev/null
@@ -0,0 +1,95 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// A test that checks that FE_Nothing does not dominate the
+// neighboring RTNodal of degree k=1,2.
+
+#include "../tests.h"
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_refinement.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/hp/fe_values.h>
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+#include <iostream>
+#include <vector>
+
+
+const unsigned int dim=2;
+
+void print_dofs (const hp::DoFHandler<2>::active_cell_iterator &cell)
+{
+  deallog << "DoFs on cell=" << cell << ": ";
+
+  std::vector<types::global_dof_index> dof_indices (cell->get_fe().dofs_per_cell);
+  cell->get_dof_indices (dof_indices);
+  for (unsigned int i=0; i<dof_indices.size(); ++i)
+    deallog << dof_indices[i] << ' ';
+  deallog << std::endl;
+}
+
+
+int main()
+{
+  initlog();
+
+  Triangulation<dim>   triangulation;
+  std::vector<unsigned int> subdivisions(dim, 1U);
+  subdivisions[0] = 2;
+  GridGenerator::subdivided_hyper_rectangle (triangulation,
+                                             subdivisions,
+                                             Point<dim>(0,0),
+                                             Point<dim>(2,1));
+  for (unsigned int i=0; i<2; ++i)
+    {
+      hp::FECollection<dim> fe_collection;
+      fe_collection.push_back(FESystem<dim>(FE_Nothing<dim>(dim),1,
+                                            FE_Nothing<dim>(),1,
+                                            FE_Q<dim>(i+1),1));
+      fe_collection.push_back(FESystem<dim>(FE_RaviartThomasNodal<dim>(i),1,
+                                            FE_DGQ<dim>(i),1,
+                                            FE_Q<dim>(i+2),1));
+
+      hp::DoFHandler<dim> dof_handler(triangulation);
+
+      dof_handler.begin_active()->set_active_fe_index(1);
+
+      dof_handler.distribute_dofs(fe_collection);
+
+      deallog << "RTNodal of degree " << i << std::endl;
+      print_dofs (dof_handler.begin_active());
+      print_dofs (++dof_handler.begin_active());
+
+      ConstraintMatrix constraints;
+      constraints.clear();
+      dealii::DoFTools::make_hanging_node_constraints  (dof_handler, constraints);
+      constraints.close ();
+
+      constraints.print(deallog.get_file_stream());
+    }
+}
+
+
+
diff --git a/tests/hp/face_domination_04.output b/tests/hp/face_domination_04.output
new file mode 100644 (file)
index 0000000..2c28d41
--- /dev/null
@@ -0,0 +1,18 @@
+
+DEAL::RTNodal of degree 0
+DEAL::DoFs on cell=0.0: 0 12 1 14 2 3 4 5 6 7 8 9 10 11 
+DEAL::DoFs on cell=0.1: 12 13 14 15 
+    4 = 0
+    5 12:  0.500000
+    5 14:  0.500000
+DEAL::RTNodal of degree 1
+DEAL::DoFs on cell=0.0: 0 30 1 32 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 
+DEAL::DoFs on cell=0.1: 30 31 32 33 34 35 36 37 38 
+    6 = 0
+    7 = 0
+    8 30:  0.323607
+    8 32:  -0.123607
+    8 34:  0.800000
+    9 30:  -0.123607
+    9 32:  0.323607
+    9 34:  0.800000

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.