]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Non local dofs implementation.
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 9 Jul 2020 15:28:27 +0000 (17:28 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Thu, 30 Jul 2020 05:26:54 +0000 (07:26 +0200)
include/deal.II/dofs/dof_accessor.h
include/deal.II/dofs/dof_accessor.templates.h
include/deal.II/dofs/dof_levels.h
include/deal.II/fe/fe_base.h
include/deal.II/fe/fe_tools.templates.h
source/dofs/dof_accessor.cc
source/fe/fe_data.cc

index 997e5befff881a923fb38ecb6a41c51766050a23..d258c0289020fa5298302e531573615b847f4dad 100644 (file)
@@ -1939,6 +1939,14 @@ public:
   void
   set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
 
+  /**
+   * Set the DoF indices of this cell to the given values. This function
+   * bypasses the DoF cache, if one exists for the given DoF handler class.
+   */
+  void
+  set_non_local_dof_indices(
+    const std::vector<types::global_dof_index> &non_local_dof_indices);
+
   /**
    * Set the Level DoF indices of this cell to the given values.
    */
index c97f88f2479f8b03c1f50c17a88b9c784bda2113..9f2d6dd98cc6e1c9a9c46a5c7eceacf35cac9efe 100644 (file)
@@ -925,7 +925,8 @@ namespace internal
       {
         AssertDimension(dof_indices.size(), n_dof_indices(accessor, fe_index));
 
-        const auto &fe = accessor.get_fe(fe_index);
+        const auto &fe                      = accessor.get_fe(fe_index);
+        const auto  non_local_dofs_per_cell = fe.n_non_local_dofs();
 
         unsigned int index = 0;
 
@@ -983,6 +984,12 @@ namespace internal
                                    dof_indices,
                                    fe_index);
 
+        // 5) non_local dofs
+        for (unsigned int d = 0; d < non_local_dofs_per_cell; ++d, ++index)
+          {
+            // TODO: Check if we need DoFOperation::process_dof also here!
+          }
+
         AssertDimension(dof_indices.size(), index);
       }
 
@@ -1682,11 +1689,51 @@ DoFAccessor<structdim, dim, spacedim, level_dof_access>::get_dof_indices(
            "been initialized, i.e., it doesn't appear that DoF indices "
            "have been distributed on it."));
 
+  switch (structdim)
+    {
+      case 1:
+        Assert(
+          dof_indices.size() ==
+            (this->n_vertices() *
+               this->dof_handler->get_fe(fe_index).n_dofs_per_vertex() +
+             this->dof_handler->get_fe(fe_index).n_dofs_per_line()) +
+              this->dof_handler->get_fe(fe_index).n_non_local_dofs_per_cell(),
+          ExcVectorDoesNotMatch());
+        break;
+      case 2:
+        Assert(
+          dof_indices.size() ==
+            (this->n_vertices() *
+               this->dof_handler->get_fe(fe_index).n_dofs_per_vertex() +
+             this->n_lines() *
+               this->dof_handler->get_fe(fe_index).n_dofs_per_line() +
+             this->dof_handler->get_fe(fe_index).n_dofs_per_quad()) +
+              this->dof_handler->get_fe(fe_index).n_non_local_dofs_per_cell(),
+          ExcVectorDoesNotMatch());
+        break;
+      case 3:
+        Assert(
+          dof_indices.size() ==
+            (this->n_vertices() *
+               this->dof_handler->get_fe(fe_index).n_dofs_per_vertex() +
+             this->n_lines() *
+               this->dof_handler->get_fe(fe_index).n_dofs_per_line() +
+             this->n_faces() *
+               this->dof_handler->get_fe(fe_index).n_dofs_per_quad() +
+             this->dof_handler->get_fe(fe_index).n_dofs_per_hex()) +
+              this->dof_handler->get_fe(fe_index).n_non_local_dofs_per_cell(),
+          ExcVectorDoesNotMatch());
+        break;
+      default:
+        Assert(false, ExcNotImplemented());
+    }
+
+
   // this function really only makes sense if either a) there are degrees of
   // freedom defined on the present object, or b) the object is non-active
-  // objects but all degrees of freedom are located on vertices, since
-  // otherwise there are degrees of freedom on sub-objects which are not
-  // allocated for this non-active thing
+  // objects but all degrees of freedom are located on vertices, since otherwise
+  // there are degrees of freedom on sub-objects which are not allocated for
+  // this non-active thing
   Assert(this->fe_index_is_active(fe_index) ||
            (this->dof_handler->get_fe(fe_index).n_dofs_per_cell() ==
             this->n_vertices() *
@@ -2205,6 +2252,126 @@ namespace internal
 
 
 
+      /**
+       * Implement setting dof indices on a cell. TO CHECK ZHAOWEI + LH
+       */
+      template <int dim, int spacedim, bool level_dof_access>
+      static void
+      set_dof_indices(
+        const DoFCellAccessor<dim, spacedim, level_dof_access> &accessor,
+        const std::vector<types::global_dof_index> &local_dof_indices)
+      {
+        Assert(accessor.has_children() == false, ExcInternalError());
+
+        const unsigned int dofs_per_vertex =
+                             accessor.get_fe().n_dofs_per_vertex(),
+                           dofs_per_line = accessor.get_fe().n_dofs_per_line(),
+                           dofs_per_quad = accessor.get_fe().n_dofs_per_quad(),
+                           dofs_per_hex  = accessor.get_fe().n_dofs_per_hex();
+
+        Assert(local_dof_indices.size() == accessor.get_fe().dofs_per_cell,
+               ExcInternalError());
+
+        unsigned int index = 0;
+
+        for (unsigned int vertex = 0;
+             vertex < GeometryInfo<dim>::vertices_per_cell;
+             ++vertex)
+          for (unsigned int d = 0; d < dofs_per_vertex; ++d, ++index)
+            accessor.set_vertex_dof_index(vertex,
+                                          d,
+                                          local_dof_indices[index],
+                                          accessor.active_fe_index());
+        // now copy dof numbers into the line. for lines in 3d with the
+        // wrong orientation, we have already made sure that we're ok
+        // by picking the correct vertices (this happens automatically
+        // in the vertex() function). however, if the line is in wrong
+        // orientation, we look at it in flipped orientation and we
+        // will have to adjust the shape function indices that we see
+        // to correspond to the correct (cell-local) ordering.
+        //
+        // of course, if dim<3, then there is nothing to adjust
+        for (unsigned int line = 0; line < GeometryInfo<dim>::lines_per_cell;
+             ++line)
+          for (unsigned int d = 0; d < dofs_per_line; ++d, ++index)
+            accessor.line(line)->set_dof_index(
+              dim < 3 ?
+                d :
+                accessor.get_fe().adjust_line_dof_index_for_line_orientation(
+                  d, accessor.line_orientation(line)),
+              local_dof_indices[index],
+              accessor.active_fe_index());
+        // now copy dof numbers into the face. for faces in 3d with the
+        // wrong orientation, we have already made sure that we're ok
+        // by picking the correct lines and vertices (this happens
+        // automatically in the line() and vertex()
+        // functions). however, if the face is in wrong orientation,
+        // we look at it in flipped orientation and we will have to
+        // adjust the shape function indices that we see to correspond
+        // to the correct (cell-local) ordering. The same applies, if
+        // the face_rotation or face_orientation is non-standard
+        //
+        // again, if dim<3, then there is nothing to adjust
+        for (unsigned int quad = 0; quad < GeometryInfo<dim>::quads_per_cell;
+             ++quad)
+          for (unsigned int d = 0; d < dofs_per_quad; ++d, ++index)
+            accessor.quad(quad)->set_dof_index(
+              dim < 3 ?
+                d :
+                accessor.get_fe().adjust_quad_dof_index_for_face_orientation(
+                  d,
+                  accessor.face_orientation(quad),
+                  accessor.face_flip(quad),
+                  accessor.face_rotation(quad)),
+              local_dof_indices[index],
+              accessor.active_fe_index());
+        for (unsigned int d = 0; d < dofs_per_hex; ++d, ++index)
+          accessor.set_dof_index(d,
+                                 local_dof_indices[index],
+                                 accessor.active_fe_index());
+        Assert(index == accessor.get_fe().dofs_per_cell, ExcInternalError());
+      }
+
+
+
+      /**
+       * Implement setting non-local dof indices on a cell.
+       */
+      template <int dim, int spacedim, bool level_dof_access>
+      static void
+      set_non_local_dof_indices(
+        const DoFCellAccessor<dim, spacedim, level_dof_access> &accessor,
+        const std::vector<types::global_dof_index> &local_non_local_dof_indices)
+      {
+        Assert(accessor.has_children() == false, ExcInternalError());
+
+        const unsigned int dofs_per_cell = accessor.get_fe().dofs_per_cell;
+
+        types::global_dof_index *next_dof_index =
+          const_cast<types::global_dof_index *>(
+            dealii::internal::DoFAccessorImplementation::Implementation::
+              get_cache_ptr(accessor.dof_handler,
+                            accessor.present_level,
+                            accessor.present_index,
+                            dofs_per_cell));
+
+        const unsigned int non_local_dofs =
+                             accessor.get_fe().n_non_local_dofs_per_cell(),
+                           n_dofs = accessor.get_fe().n_dofs_per_cell();
+
+        unsigned int index = 0;
+
+        for (unsigned int d = 0; d < n_dofs; ++d, ++next_dof_index)
+          if (d >= n_dofs - non_local_dofs)
+            {
+              *next_dof_index = local_non_local_dof_indices[index];
+              ++index;
+            }
+        Assert(index == accessor.get_fe().n_non_local_dofs_per_cell(),
+               ExcInternalError());
+      }
+
+
       /**
        * Do what the active_fe_index function in the parent class is supposed to
        * do.
index 4512e520e19996ac44cc6ce6212015018ef67f84..8170dc23fd58a1638f451fdd7592db9fe008cdef 100644 (file)
@@ -80,6 +80,12 @@ namespace internal
        */
       DoFObjects<dim> dof_object;
 
+      /**
+       * The offsets for each cell of the cache that holds all DoF indices.
+       * Only used for hp elements.
+       */
+      std::vector<unsigned int> cell_cache_offsets;
+
       /**
        * Return a pointer to the beginning of the DoF indices cache for a
        * given cell.
index 3b69a65a82913717375454f3130a35f3c02ff435..d89a618c895b8ae7c56e447be67acb14051c92c8 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2000 - 2018 by the deal.II authors
+// Copyright (C) 2000 - 2020 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -249,6 +249,23 @@ public:
    */
   const unsigned int dofs_per_hex;
 
+  /**
+   * Number of non-zero degrees of freedom (DoFs) on the cell that,
+   * however, can not be associated with the specific geometric object (such as
+   * face or edge). In these cases, DoFs are typically globally distributed over
+   * (patches of) the Triangulation. Typical examples are:
+   *  - <a
+   * href="https://www.sciencedirect.com/science/article/pii/S0045782504005171">Isogeometric
+   * Analysis</a>
+   *  - <a
+   * href="https://onlinelibrary.wiley.com/doi/abs/10.1002/%28SICI%291097-0207%2819990910%2946%3A1%3C131%3A%3AAID-NME726%3E3.0.CO%3B2-J">Extended
+   * finite element method</a>
+   *  - <a
+   * href="http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.21.4515&rep=rep1&type=pdf">Catmull-Clark's
+   * subdivision surfaces</a>
+   */
+  const unsigned int non_local_dofs_per_cell;
+
   /**
    * First index of dof on a line.
    */
@@ -321,12 +338,16 @@ public:
    * to geometrical objects.
    *
    * @param[in] dofs_per_object A vector that describes the number of degrees
-   * of freedom on geometrical objects for each dimension. This vector must
-   * have size dim+1, and entry 0 describes the number of degrees of freedom
-   * per vertex, entry 1 the number of degrees of freedom per line, etc. As an
-   * example, for the common $Q_1$ Lagrange element in 2d, this vector would
-   * have elements <code>(1,0,0)</code>. On the other hand, for a $Q_3$
-   * element in 3d, it would have entries <code>(1,2,4,8)</code>.
+   * of freedom on geometrical objects for each dimension.
+   * The first dim+1 elements of the vector describe the number
+   * of degrees of freedom per geometrical object: entry 0 describes the number
+   * of degrees of freedom per vertex, entry 1 the number of degrees of freedom
+   * per line, etc. If the vector is of size dim+2, the last entry specifices
+   * the number of non-local degrees of freedom for this element. In case the
+   * vector size is dim+1, the number of non-local degrees of freedom is set to
+   * zero. As an example, for the common $Q_1$ Lagrange element in 2d, this
+   * vector would have elements <code>(1,0,0)</code>. On the other hand, for a
+   * $Q_3$ element in 3d, it would have entries <code>(1,2,4,8)</code>.
    *
    * @param[in] n_components Number of vector components of the element.
    *
@@ -415,6 +436,24 @@ public:
   unsigned int
   n_dofs_per_cell() const;
 
+  /**
+   * Return the number of non-zero degrees of freedom (DoFs) on the cell that,
+   * however, can not be associated with the specific geometric object (such as
+   * face or edge). In these cases, DoFs are typically globally distributed over
+   * (patches of) the Triangulation. Typical examples are:
+   *  - <a
+   * href="https://www.sciencedirect.com/science/article/pii/S0045782504005171">Isogeometric
+   * Analysis</a>
+   *  - <a
+   * href="https://onlinelibrary.wiley.com/doi/abs/10.1002/%28SICI%291097-0207%2819990910%2946%3A1%3C131%3A%3AAID-NME726%3E3.0.CO%3B2-J">Extended
+   * finite element method</a>
+   *  - <a
+   * href="http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.21.4515&rep=rep1&type=pdf">Catmull-Clark's
+   * subdivision surfaces</a>
+   */
+  unsigned int
+  n_non_local_dofs_per_cell() const;
+
   /**
    * Return the number of degrees per structdim-dimensional object. For
    * structdim==0, the function therefore returns dofs_per_vertex, for
@@ -621,6 +660,15 @@ FiniteElementData<dim>::n_dofs_per_cell() const
 
 
 
+template <int dim>
+inline unsigned int
+FiniteElementData<dim>::n_non_local_dofs_per_cell() const
+{
+  return non_local_dofs_per_cell;
+}
+
+
+
 template <int dim>
 template <int structdim>
 inline unsigned int
index 34c7d92789867daab53bbec0e2aa2ebc338bbc39..046cd7d6c5755ec228ee80f25ebe1efe8e5ad7e4 100644 (file)
@@ -85,10 +85,11 @@ namespace FETools
     {
       AssertDimension(fes.size(), multiplicities.size());
 
-      unsigned int multiplied_dofs_per_vertex = 0;
-      unsigned int multiplied_dofs_per_line   = 0;
-      unsigned int multiplied_dofs_per_quad   = 0;
-      unsigned int multiplied_dofs_per_hex    = 0;
+      unsigned int multiplied_dofs_per_vertex         = 0;
+      unsigned int multiplied_dofs_per_line           = 0;
+      unsigned int multiplied_dofs_per_quad           = 0;
+      unsigned int multiplied_dofs_per_hex            = 0;
+      unsigned int multiplied_non_local_dofs_per_cell = 0;
 
       unsigned int multiplied_n_components = 0;
 
@@ -114,6 +115,8 @@ namespace FETools
               fes[i]->n_dofs_per_quad() * multiplicities[i];
             multiplied_dofs_per_hex +=
               fes[i]->n_dofs_per_hex() * multiplicities[i];
+            multiplied_non_local_dofs_per_cell +=
+              fes[i]->n_non_local_dofs_per_cell() * multiplicities[i];
 
             multiplied_n_components +=
               fes[i]->n_components() * multiplicities[i];
@@ -152,6 +155,8 @@ namespace FETools
         dpo.push_back(multiplied_dofs_per_quad);
       if (dim > 2)
         dpo.push_back(multiplied_dofs_per_hex);
+      if (multiplied_non_local_dofs_per_cell > 0)
+        dpo.push_back(multiplied_non_local_dofs_per_cell);
 
       BlockIndices block_indices(0, 0);
 
@@ -319,7 +324,36 @@ namespace FETools
                     fes[base]->restriction_is_additive(index_in_base);
                 }
         }
+      unsigned int total_index_offset = total_index;
+      // 5. Non_local
+      unsigned int n_fe = 0;
+      for (unsigned int base = 0; base < fes.size(); ++base)
+        {
+          if (fes[base] != NULL && fes[base]->n_non_local_dofs_per_cell() > 0)
+            {
+              ++n_fe;
+            }
+        }
+      for (unsigned int base = 0; base < n_fe; ++base)
+        for (unsigned int m = 0; m < multiplicities[base]; ++m)
+          for (unsigned int local_index = 0;
+               local_index < fes[base]->n_non_local_dofs_per_cell();
+               ++local_index)
+            {
+              total_index = total_index_offset +
+                            local_index * (n_fe * multiplicities[base]) +
+                            base * multiplicities[base] + m;
+
+              const unsigned int index_in_base =
+                (local_index + (fes[base]->n_dofs_per_cell() -
+                                fes[base]->n_non_local_dofs_per_cell()));
+              Assert(index_in_base < fes[base]->n_dofs_per_cell(),
+                     ExcInternalError());
+              retval[total_index] =
+                fes[base]->restriction_is_additive(index_in_base);
+            }
 
+      // TODO[LH]: Make sure this Assert here is correct
       Assert(total_index == n_shape_functions, ExcInternalError());
 
       return retval;
@@ -378,8 +412,8 @@ namespace FETools
     {
       AssertDimension(fes.size(), multiplicities.size());
 
-      // first count the number of dofs and components that will emerge from the
-      // given FEs
+      // first count the number of dofs and components that will emerge
+      // from the given FEs
       unsigned int n_shape_functions = 0;
       for (unsigned int i = 0; i < fes.size(); ++i)
         if (multiplicities[i] > 0) // needed because fe might be nullptr
@@ -413,15 +447,15 @@ namespace FETools
                                             std::vector<bool>(n_components,
                                                               false));
 
-      // finally go through all the shape functions of the base elements, and
-      // copy their flags. this somehow copies the code in build_cell_table,
-      // which is not nice as it uses too much implicit knowledge about the
-      // layout of the individual bases in the composed FE, but there seems no
-      // way around...
+      // finally go through all the shape functions of the base elements,
+      // and copy their flags. this somehow copies the code in
+      // build_cell_table, which is not nice as it uses too much implicit
+      // knowledge about the layout of the individual bases in the
+      // composed FE, but there seems no way around...
       //
-      // for each shape function, copy the non-zero flags from the base element
-      // to this one, taking into account multiplicities, multiple components in
-      // base elements, and other complications
+      // for each shape function, copy the non-zero flags from the base
+      // element to this one, taking into account multiplicities, multiple
+      // components in base elements, and other complications
       unsigned int total_index = 0;
       for (const unsigned int vertex_number :
            ReferenceCell::internal::Info::get_cell(
@@ -561,7 +595,54 @@ namespace FETools
                     }
                 }
         }
+      unsigned int total_index_offset = total_index;
+      // 5. non_local
+      unsigned int comp_start = 0;
+      unsigned int n_fe       = 0;
+      for (unsigned int base = 0; base < fes.size(); ++base)
+        {
+          if (fes[base] != NULL && fes[base]->n_non_local_dofs_per_cell() > 0)
+            {
+              ++n_fe;
+            }
+        }
+      for (unsigned int base = 0; base < n_fe; ++base)
+        {
+          for (unsigned int m = 0; m < multiplicities[base];
+               ++m, comp_start += fes[base]->n_components() * do_tensor_product)
+            {
+              for (unsigned int non_local_index = 0;
+                   non_local_index < fes[base]->n_non_local_dofs_per_cell();
+                   ++non_local_index)
+                {
+                  total_index =
+                    total_index_offset +
+                    non_local_index * (n_fe * multiplicities[base]) +
+                    base * multiplicities[base] + m;
+
+                  const unsigned int index_in_base =
+                    (non_local_index +
+                     (fes[base]->n_dofs_per_cell() -
+                      fes[base]->n_non_local_dofs_per_cell()));
+
+                  Assert(comp_start + fes[base]->n_components() <=
+                           retval[total_index].size(),
+                         ExcInternalError());
+
+                  for (unsigned int c = 0; c < fes[base]->n_components(); ++c)
+                    {
+                      Assert(c < fes[base]
+                                   ->get_nonzero_components(index_in_base)
+                                   .size(),
+                             ExcInternalError());
+                      retval[total_index][comp_start + c] =
+                        fes[base]->get_nonzero_components(index_in_base)[c];
+                    }
+                }
+            }
+        }
 
+      // TODO[LH]: Make sure this assert here is correct
       Assert(total_index == n_shape_functions, ExcInternalError());
 
       // now copy the vector<vector<bool> > into a vector<ComponentMask>.
@@ -843,6 +924,57 @@ namespace FETools
                       non_primitive_index;
                 }
         }
+      unsigned int total_index_offset = total_index;
+      // 5. non-local
+      for (unsigned int base = 0; base < fe.n_base_elements(); ++base)
+        {
+          unsigned int comp_start = 0;
+          for (unsigned int m = 0; m < fe.element_multiplicity(base);
+               ++m,
+                            comp_start += fe.base_element(base).n_components() *
+                                          do_tensor_product)
+            {
+              for (unsigned int non_local_index = 0;
+                   non_local_index <
+                   fe.base_element(base).n_non_local_dofs_per_cell();
+                   ++non_local_index)
+                {
+                  total_index =
+                    total_index_offset +
+                    non_local_index *
+                      (fe.n_base_elements() * fe.element_multiplicity(base)) +
+                    base * fe.element_multiplicity(base) + m;
+
+                  const unsigned int index_in_base =
+                    (non_local_index +
+                     (fe.base_element(base).n_dofs_per_cell() -
+                      fe.base_element(base).n_non_local_dofs_per_cell()));
+
+                  system_to_base_table[total_index] =
+                    std::make_pair(std::make_pair(base, m), index_in_base);
+
+                  if (fe.base_element(base).is_primitive(index_in_base))
+                    {
+                      const unsigned int comp_in_base =
+                        fe.base_element(base)
+                          .system_to_component_index(index_in_base)
+                          .first;
+                      const unsigned int comp = comp_start + comp_in_base;
+                      const unsigned int index_in_comp =
+                        fe.base_element(base)
+                          .system_to_component_index(index_in_base)
+                          .second;
+                      system_to_component_table[total_index] =
+                        std::make_pair(comp, index_in_comp);
+                    }
+                  else
+                    {
+                      system_to_component_table[total_index] =
+                        non_primitive_index;
+                    }
+                }
+            }
+        }
     }
 
 
@@ -885,12 +1017,12 @@ namespace FETools
                 {
                   // get (cell) index of this shape function inside the base
                   // element to see whether the shape function is primitive
-                  // (assume that all shape functions on vertices share the same
-                  // primitivity property; assume likewise for all shape
+                  // (assume that all shape functions on vertices share the
+                  // same primitivity property; assume likewise for all shape
                   // functions located on lines, quads, etc. this way, we can
-                  // ask for primitivity of only _one_ shape function, which is
-                  // taken as representative for all others located on the same
-                  // type of object):
+                  // ask for primitivity of only _one_ shape function, which
+                  // is taken as representative for all others located on the
+                  // same type of object):
                   const unsigned int index_in_base =
                     (fe.base_element(base).n_dofs_per_vertex() * vertex_number +
                      local_index);
index 579536a127e9a560384aaf6b86dfbc73754eea00..4f649e9906d23145c2e009a6a12f1b3740be638c 100644 (file)
@@ -138,6 +138,23 @@ DoFCellAccessor<dim, spacedim, lda>::set_dof_indices(
 
 
 
+template <int dim, int spacedim, bool lda>
+void
+DoFCellAccessor<dim, spacedim, lda>::set_non_local_dof_indices(
+  const std::vector<types::global_dof_index> &local_non_local_dof_indices)
+{
+  Assert(static_cast<unsigned int>(this->present_level) <
+           this->dof_handler->object_dof_indices.size(),
+         ExcMessage("DoFHandler not initialized"));
+
+  Assert(this->dof_handler != nullptr, typename BaseClass::ExcInvalidObject());
+
+  internal::DoFCellAccessorImplementation::Implementation::
+    set_non_local_dof_indices(*this, local_non_local_dof_indices);
+}
+
+
+
 template <int dim, int spacedim, bool lda>
 TriaIterator<DoFCellAccessor<dim, spacedim, lda>>
 DoFCellAccessor<dim, spacedim, lda>::neighbor_child_on_subface(
index f9fff482ff17ff37efd55e5a4f244d07f196ff22..81ebd0bb65b4b22a9a7b4faf41c0edb396c98a70 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2001 - 2018 by the deal.II authors
+// Copyright (C) 2001 - 2020 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -51,6 +51,8 @@ FiniteElementData<dim>::FiniteElementData(
   , dofs_per_line(dofs_per_object[1])
   , dofs_per_quad(dim > 1 ? dofs_per_object[2] : 0)
   , dofs_per_hex(dim > 2 ? dofs_per_object[3] : 0)
+  , non_local_dofs_per_cell(
+      dofs_per_object.size() == dim + 2 ? dofs_per_object[dim + 1] : 0)
   , first_line_index(
       ReferenceCell::internal::Info::get_cell(cell_type).n_vertices() *
       dofs_per_vertex)
@@ -94,7 +96,7 @@ FiniteElementData<dim>::FiniteElementData(
             ReferenceCell::internal::Info::get_cell(cell_type).n_faces() :
             0)) *
         dofs_per_quad +
-      (dim == 3 ? 1 : 0) * dofs_per_hex)
+      (dim == 3 ? 1 : 0) * dofs_per_hex + non_local_dofs_per_cell)
   , components(n_components)
   , degree(degree)
   , conforming_space(conformity)
@@ -102,8 +104,10 @@ FiniteElementData<dim>::FiniteElementData(
                          BlockIndices(1, dofs_per_cell) :
                          block_indices)
 {
-  Assert(dofs_per_object.size() == dim + 1,
-         ExcDimensionMismatch(dofs_per_object.size() - 1, dim));
+  Assert(dofs_per_object.size() == dim + 1 || dofs_per_object.size() == dim + 2,
+         ExcMessage("dofs_per_object should have size of either " +
+                    std::to_string(dim + 1) + " or " +
+                    std::to_string(dim + 2)));
 }
 
 
@@ -116,7 +120,8 @@ FiniteElementData<dim>::operator==(const FiniteElementData<dim> &f) const
           (dofs_per_line == f.dofs_per_line) &&
           (dofs_per_quad == f.dofs_per_quad) &&
           (dofs_per_hex == f.dofs_per_hex) && (components == f.components) &&
-          (degree == f.degree) && (conforming_space == f.conforming_space));
+          (degree == f.degree) && (conforming_space == f.conforming_space) &&
+          (non_local_dofs_per_cell == f.non_local_dofs_per_cell));
 }
 
 

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.