Non local dofs within FE.

diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h
index 08a9006f025059810488b9a655e70afd2b856fcb..bef69bba72726c83455d814831d3addf23ca48b1 100644 (file)
--- a/include/deal.II/dofs/dof_accessor.templates.h
+++ b/include/deal.II/dofs/dof_accessor.templates.h
@@ -822,11 +822,12 @@ namespace internal
           {
             const auto &fe = accessor.get_fe(fe_index_);
 
-            const unsigned int                          //
-              dofs_per_vertex = fe.n_dofs_per_vertex(), //
-              dofs_per_line   = fe.n_dofs_per_line(),   //
-              dofs_per_quad   = fe.n_dofs_per_quad(),   //
-              dofs_per_hex    = fe.n_dofs_per_hex();    //
+            const unsigned int                                  //
+              dofs_per_vertex = fe.n_dofs_per_vertex(),         //
+              dofs_per_line   = fe.n_dofs_per_line(),           //
+              dofs_per_quad   = fe.n_dofs_per_quad(),           //
+              dofs_per_hex    = fe.n_dofs_per_hex(),            //
+              non_local_dofs  = fe.n_non_local_dofs_per_cell(); //
 
             const unsigned int inner_dofs =
               structdim == 1 ? dofs_per_line :
@@ -860,6 +861,10 @@ namespace internal
             // 4) INNER dofs
             index += inner_dofs;
 
+            // 5) non local dofs
+            if (dim == structdim)
+              index += non_local_dofs;
+
             return index;
           }
         else
@@ -874,6 +879,9 @@ namespace internal
 
             const auto diff = [](const auto &p) { return p.second - p.first; };
 
+            const unsigned int non_local_dofs =
+              accessor.get_fe(fe_index).n_non_local_dofs_per_cell();
+
             // 1) VERTEX dofs
             for (const auto vertex : accessor.vertex_indices())
               index +=
@@ -898,6 +906,10 @@ namespace internal
             // 4) INNER dofs
             index += diff(process_object_range(accessor, fe_index));
 
+            // 5) non local dofs
+            if (dim == structdim)
+              index += non_local_dofs;
+
             return index;
           }
       }
@@ -926,7 +938,7 @@ namespace internal
         AssertDimension(dof_indices.size(), n_dof_indices(accessor, fe_index));
 
         const auto &fe                      = accessor.get_fe(fe_index);
-        const auto  non_local_dofs_per_cell = fe.n_non_local_dofs();
+        const auto  non_local_dofs_per_cell = fe.n_non_local_dofs_per_cell();
 
         unsigned int index = 0;
 

diff --git a/include/deal.II/fe/fe.h b/include/deal.II/fe/fe.h
index 5ccb1589f9396320056d555f28c1136c1ef3c6df..7882b2d05e03f6720757d2e135ab8d5f1284fbbb 100644 (file)
--- a/include/deal.II/fe/fe.h
+++ b/include/deal.II/fe/fe.h
@@ -44,6 +44,9 @@ template <int dim, int spacedim>
 class FESubfaceValues;
 template <int dim, int spacedim>
 class FESystem;
+template <int dim, int spacedim, bool>
+class DoFCellAccessor;
+
 
 /**
  * This is the base class for finite elements in arbitrary dimensions. It
@@ -2267,10 +2270,50 @@ public:
    * @{
    */
   /**
-   * Return an object that knows how to handle non local dof indices.
+   * Return the non local dof indices associated to the current cell, for
+   * active cell accessors.
+   */
+  virtual std::vector<types::global_dof_index>
+  get_non_local_dof_indices(
+    const DoFCellAccessor<dim, spacedim, false> &accessor) const;
+
+  /**
+   * Return the non local dof indices associated to the current cell, for
+   * level cell accessors.
+   */
+  virtual std::vector<types::global_dof_index>
+  get_non_local_dof_indices(
+    const DoFCellAccessor<dim, spacedim, true> &accessor) const;
+
+  /**
+   * Return the global number of non local dof indices that are required in
+   * addition to the local ones.
+   */
+  virtual types::global_dof_index
+  n_global_non_local_dofs() const;
+
+  /**
+   * Return an identification string that uniquely identifies the non local
+   * behaviour of the finite element space.
+   *
+   * For non local finite element spaces, n_non_local_dofs_per_cell() and
+   * n_global_non_local_dofs() may return a non zero number, meaning that there
+   * are degrees of freedom that are not associated to vertices, edges, faces,
+   * or cells (for example, they may be associated to patches of cells, or be
+   * global non zero basis functions), that are non zero on certain cells.
+   *
+   * In these cases, one usually uses the hp support of the library, and defines
+   * an FECollection where each FiniteElement of the collection has the same non
+   * local behaviour. This id is checked when calling
+   * DoFHandler::distribute_dofs() with an FECollection as argument, and an
+   * assertion is thrown if the ids do not coincide for all FiniteElement
+   * objects of the collection.
+   *
+   * By default, FiniteElement spaces are local, and this function returns the
+   * string "Local FiniteElement space".
    */
-  virtual std::shared_ptr<const NonLocalDoFHandler<dim, spacedim>>
-  get_non_local_dof_handler() const;
+  virtual std::string
+  get_non_local_id() const;
   //@}
 
   /**
@@ -3314,14 +3357,52 @@ FiniteElement<dim, spacedim>::get_associated_geometry_primitive(
 
 
 
+// Non local dofs support.
+template <int dim, int spacedim>
+inline std::vector<types::global_dof_index>
+FiniteElement<dim, spacedim>::get_non_local_dof_indices(
+  const DoFCellAccessor<dim, spacedim, true> &) const
+{
+  Assert(this->n_non_local_dofs_per_cell() == 0 &&
+           n_global_non_local_dofs() == 0,
+         ExcPureFunctionCalled());
+  return std::vector<types::global_dof_index>();
+}
+
+
+
+template <int dim, int spacedim>
+inline std::vector<types::global_dof_index>
+FiniteElement<dim, spacedim>::get_non_local_dof_indices(
+  const DoFCellAccessor<dim, spacedim, false> &) const
+{
+  Assert(this->n_non_local_dofs_per_cell() == 0 &&
+           n_global_non_local_dofs() == 0,
+         ExcPureFunctionCalled());
+  return std::vector<types::global_dof_index>();
+}
+
+
+
 template <int dim, int spacedim>
-inline std::shared_ptr<const NonLocalDoFHandler<dim, spacedim>>
-FiniteElement<dim, spacedim>::get_non_local_dof_handler() const
+inline types::global_dof_index
+FiniteElement<dim, spacedim>::n_global_non_local_dofs() const
 {
-  return std::make_shared<const NonLocalDoFHandler<dim, spacedim>>();
+  return 0;
 }
 
 
+
+template <int dim, int spacedim>
+inline std::string
+FiniteElement<dim, spacedim>::get_non_local_id() const
+{
+  Assert(this->n_non_local_dofs_per_cell() == 0 &&
+           n_global_non_local_dofs() == 0,
+         ExcPureFunctionCalled());
+  return "Local FiniteElement space";
+}
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif

diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc
index 1485977780ed7d3950e3341510e0356f9ab72e4e..9f5c6a07f4bc948f65e7e8679e5d11e177c247fd 100644 (file)
--- a/source/dofs/dof_handler.cc
+++ b/source/dofs/dof_handler.cc
@@ -2612,18 +2612,22 @@ DoFHandler<dim, spacedim>::distribute_non_local_dofs()
     ExcMessage(
       "Distribute active DoFs using distribute_dofs() before calling distribute_non_local_dofs()."));
 
-  auto non_local_dh = begin()->get_fe().get_non_local_dof_handler();
+  const auto non_local_id = begin()->get_fe().get_non_local_id();
+  const auto n_global_non_local_dofs =
+    begin()->get_fe().n_global_non_local_dofs();
+
   for (auto cell : active_cell_iterators())
     {
-      // Make sure all cells use the same NonLocalDoFHandler.
-      Assert(non_local_dh == cell->get_fe().get_non_local_dof_handler(),
-             ExcMessage(
-               "You are trying to use different NonLocalDoFHandler objects"));
-      cell->set_non_local_dof_indices(
-        non_local_dh->get_non_local_dof_indices(*cell));
+      const auto &fe = cell->get_fe();
+      // Make sure all cells handle non local dofs in the same way.
+      Assert(
+        non_local_id == fe.get_non_local_id(),
+        ExcMessage(
+          "You are trying to use different non local finite element spaces"));
+      AssertDimension(n_global_non_local_dofs, fe.n_global_non_local_dofs());
+      cell->set_non_local_dof_indices(fe.get_non_local_dof_indices(*cell));
     }
-  this->number_cache.n_global_dofs +=
-    non_local_dh->n_additional_non_local_dofs();
+  this->number_cache.n_global_dofs += n_global_non_local_dofs;
 }
 
 

diff --git a/tests/non_local/fe_q1_nonlocal.h b/tests/non_local/fe_q1_nonlocal.h
index 1d713624ac3841b51a98ce8534f981aa48035ead..f9a7c1c1111f661526be238e4df16ca331900854 100644 (file)
--- a/tests/non_local/fe_q1_nonlocal.h
+++ b/tests/non_local/fe_q1_nonlocal.h
@@ -21,6 +21,8 @@
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/tensor_product_polynomials.h>
 
+#include <deal.II/dofs/dof_accessor.h>
+
 #include <deal.II/fe/fe_q_base.h>
 
 #include <iostream>
@@ -36,45 +38,11 @@ get_dpo_vector()
   return dpo;
 }
 
-
-/**
- * Associates non local dofs according to the vertex index.
- */
-template <int dim>
-class NonLocalQ1DoFHandler : public NonLocalDoFHandler<dim>
-{
-public:
-  virtual std::vector<types::global_dof_index>
-  get_non_local_dof_indices(
-    const DoFCellAccessor<dim, dim, false> &accessor) const override
-  {
-    if (!tria)
-      tria = &(accessor.get_triangulation());
-    std::vector<types::global_dof_index> dofs(accessor.n_vertices());
-    for (unsigned int i = 0; i < dofs.size(); ++i)
-      dofs[i] = accessor.vertex_index(i);
-    return dofs;
-  }
-
-  virtual types::global_dof_index
-  n_additional_non_local_dofs() const override
-  {
-    Assert(tria, ExcInternalError());
-    return tria->n_vertices();
-  }
-
-private:
-  mutable SmartPointer<const Triangulation<dim>> tria = nullptr;
-};
-
-
-
 template <int dim>
 class FE_Q1_Nonlocal : public FE_Q_Base<TensorProductPolynomials<dim>, dim, dim>
 {
 public:
-  FE_Q1_Nonlocal(const std::shared_ptr<NonLocalQ1DoFHandler<dim>> &ptr =
-                   std::shared_ptr<NonLocalQ1DoFHandler<dim>>())
+  FE_Q1_Nonlocal(const Triangulation<dim> &tria)
     : FE_Q_Base<TensorProductPolynomials<dim>, dim, dim>(
         TensorProductPolynomials<dim>(
           Polynomials::generate_complete_Lagrange_basis(
@@ -84,7 +52,7 @@ public:
                                1,
                                FiniteElementData<dim>::H1),
         std::vector<bool>(1, false))
-    , non_local_dh(ptr ? ptr : std::make_shared<NonLocalQ1DoFHandler<dim>>())
+    , tria(&tria)
   {
     this->unit_support_points = QTrapez<dim>().get_points();
   }
@@ -92,7 +60,7 @@ public:
   virtual std::unique_ptr<FiniteElement<dim>>
   clone() const override
   {
-    return std::make_unique<FE_Q1_Nonlocal<dim>>(non_local_dh);
+    return std::make_unique<FE_Q1_Nonlocal<dim>>(*tria);
   }
 
   virtual std::string
@@ -101,12 +69,6 @@ public:
     return "FE_Q_Nonlocal<dim>";
   }
 
-  virtual std::shared_ptr<const NonLocalDoFHandler<dim>>
-  get_non_local_dof_handler() const override
-  {
-    return non_local_dh;
-  }
-
   virtual void
   convert_generalized_support_point_values_to_dof_values(
     const std::vector<Vector<double>> &support_point_values,
@@ -124,6 +86,31 @@ public:
       }
   }
 
+
+  virtual std::vector<types::global_dof_index>
+  get_non_local_dof_indices(
+    const DoFCellAccessor<dim, dim, false> &accessor) const override
+  {
+    std::vector<types::global_dof_index> dofs(accessor.n_vertices());
+    for (unsigned int i = 0; i < dofs.size(); ++i)
+      dofs[i] = accessor.vertex_index(i);
+    return dofs;
+  }
+
+
+  virtual types::global_dof_index
+  n_global_non_local_dofs() const override
+  {
+    Assert(tria, ExcInternalError());
+    return tria->n_vertices();
+  }
+
+  virtual std::string
+  get_non_local_id() const override
+  {
+    return "NonLocal FEQ1";
+  }
+
 private:
-  std::shared_ptr<NonLocalQ1DoFHandler<dim>> non_local_dh;
+  mutable SmartPointer<const Triangulation<dim>> tria = nullptr;
 };

diff --git a/tests/non_local/fe_q1_nonlocal_01.cc b/tests/non_local/fe_q1_nonlocal_01.cc
index 5633da7015414472dfad308411190139994b739d..56f2c73c46b5fe23f86c60fb04545975f5f3ba3f 100644 (file)
--- a/tests/non_local/fe_q1_nonlocal_01.cc
+++ b/tests/non_local/fe_q1_nonlocal_01.cc
@@ -31,7 +31,8 @@ template <int dim>
 void
 test()
 {
-  FE_Q1_Nonlocal<dim> fe;
+  Triangulation<dim>  tria;
+  FE_Q1_Nonlocal<dim> fe(tria);
   if (fe.n_dofs_per_cell() != GeometryInfo<dim>::vertices_per_cell)
     deallog << "Not OK" << std::endl;
   else

diff --git a/tests/non_local/fe_q1_nonlocal_02.cc b/tests/non_local/fe_q1_nonlocal_02.cc
index a490caa19236173d5bed95da2bc0a0c652e68501..0339d6164232af772961517626c6013eb56892d8 100644 (file)
--- a/tests/non_local/fe_q1_nonlocal_02.cc
+++ b/tests/non_local/fe_q1_nonlocal_02.cc
@@ -31,7 +31,8 @@ void
 test()
 {
   deallog << "DIM: " << dim << std::endl << std::endl;
-  FE_Q1_Nonlocal<dim> fe;
+  Triangulation<dim>  tria;
+  FE_Q1_Nonlocal<dim> fe(tria);
   // Check that the basis functions are what we expect
   const auto &quad = fe.get_unit_support_points();
 

diff --git a/tests/non_local/fe_q1_nonlocal_03.cc b/tests/non_local/fe_q1_nonlocal_03.cc
index 5746bba08866eff2aae626aa64593c48da1d6df9..8a06b270b2b1645ca835649cfe95cb5ac85a295c 100644 (file)
--- a/tests/non_local/fe_q1_nonlocal_03.cc
+++ b/tests/non_local/fe_q1_nonlocal_03.cc
@@ -31,9 +31,9 @@ template <int dim>
 void
 test()
 {
-  FE_Q1_Nonlocal<dim> fe;
-  Triangulation<dim>  tria;
+  Triangulation<dim> tria;
   GridGenerator::hyper_cube(tria);
+  FE_Q1_Nonlocal<dim> fe(tria);
   tria.refine_global(1);
   DoFHandler<dim> dh(tria);
   dh.distribute_dofs(fe);

diff --git a/tests/non_local/fe_q1_nonlocal_04.cc b/tests/non_local/fe_q1_nonlocal_04.cc
index cebe016ec82fd88f6b62e3360d5aeefd58c68810..5150249355aca70dab7549dd56c63fa64255cb76 100644 (file)
--- a/tests/non_local/fe_q1_nonlocal_04.cc
+++ b/tests/non_local/fe_q1_nonlocal_04.cc
@@ -35,9 +35,9 @@ template <int dim>
 void
 test()
 {
-  FE_Q1_Nonlocal<dim> fe;
-  Triangulation<dim>  tria;
+  Triangulation<dim> tria;
   GridGenerator::hyper_cube(tria);
+  FE_Q1_Nonlocal<dim> fe(tria);
   tria.refine_global(1);
   DoFHandler<dim> dh(tria);
   dh.distribute_dofs(fe);