]> https://gitweb.dealii.org/ - dealii.git/commitdiff
periodic_neighbor in DoFTools::make_sparsity_pattern
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 29 Mar 2016 17:26:54 +0000 (19:26 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 18 Apr 2016 17:10:19 +0000 (19:10 +0200)
include/deal.II/dofs/dof_accessor.h
source/dofs/dof_accessor.cc
source/dofs/dof_tools_sparsity.cc

index e25e74412fe06d785d6658ba04c4709dc809e094..8b257ed11e75dd22b644da5e2432d685207db05a 100644 (file)
@@ -1264,6 +1264,22 @@ public:
   TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
   neighbor (const unsigned int) const;
 
+  /**
+   * Return the @p ith periodic neighbor as a DoF cell iterator. This function
+   * is needed since the neighbor function of the base class returns a cell
+   * accessor without access to the DoF data.
+   */
+  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
+  periodic_neighbor (const unsigned int) const;
+
+  /**
+   * Return the @p ith neighbor or periodic neighbor as a DoF cell iterator.
+   * This function is needed since the neighbor function of the base class
+   * returns a cell accessor without access to the DoF data.
+   */
+  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
+  neighbor_or_periodic_neighbor (const unsigned int) const;
+
   /**
    * Return the @p ith child as a DoF cell iterator. This function is needed
    * since the child function of the base class returns a cell accessor
@@ -1291,6 +1307,16 @@ public:
   neighbor_child_on_subface (const unsigned int face_no,
                              const unsigned int subface_no) const;
 
+  /**
+   * Return the result of the @p periodic_neighbor_child_on_subface function
+   * of the base class, but convert it so that one can also access the DoF
+   * data (the function in the base class only returns an iterator with access
+   * to the triangulation data).
+   */
+  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
+  periodic_neighbor_child_on_subface (const unsigned int face_no,
+                                      const unsigned int subface_no) const;
+
   /**
    * @}
    */
index 588654d0d9a7b13f06235ebc48f66878d7ca0cae..b2751a96637c52b0e5d3c9c0e928fa35bee650f0 100644 (file)
@@ -84,6 +84,43 @@ DoFCellAccessor<DoFHandlerType,lda>::neighbor_child_on_subface (const unsigned i
 
 
 
+template <typename DoFHandlerType, bool lda>
+TriaIterator<DoFCellAccessor<DoFHandlerType,lda> >
+DoFCellAccessor<DoFHandlerType,lda>::periodic_neighbor_child_on_subface (const unsigned int face,
+    const unsigned int subface) const
+{
+  const TriaIterator<CellAccessor<dim,spacedim> > q
+    = CellAccessor<dim,spacedim>::periodic_neighbor_child_on_subface (face, subface);
+  return TriaIterator<DoFCellAccessor<DoFHandlerType,lda> > (*q, this->dof_handler);
+}
+
+
+
+template <typename DoFHandlerType, bool lda>
+inline
+TriaIterator<DoFCellAccessor<DoFHandlerType,lda> >
+DoFCellAccessor<DoFHandlerType,lda>::periodic_neighbor (const unsigned int face) const
+{
+  const TriaIterator<CellAccessor<dim,spacedim> > q
+    = CellAccessor<dim,spacedim>::periodic_neighbor (face);
+  return TriaIterator<DoFCellAccessor<DoFHandlerType,lda> > (*q, this->dof_handler);
+}
+
+
+
+template <typename DoFHandlerType, bool lda>
+inline
+TriaIterator<DoFCellAccessor<DoFHandlerType,lda> >
+DoFCellAccessor<DoFHandlerType,lda>::neighbor_or_periodic_neighbor (const unsigned int face) const
+{
+  const TriaIterator<CellAccessor<dim,spacedim> > q
+    = CellAccessor<dim,spacedim>::neighbor_or_periodic_neighbor (face);
+  return TriaIterator<DoFCellAccessor<DoFHandlerType,lda> > (*q, this->dof_handler);
+}
+
+
+
+
 // --------------------------------------------------------------------------
 // explicit instantiations
 #include "dof_accessor.inst"
index 5a98b4e698b58a832ec0adf228c46b9d95ad2ae9..7d9db73589f05e84e9c94c382e41e8374db65beb 100644 (file)
@@ -553,13 +553,15 @@ namespace DoFTools
                ++face)
             {
               typename DoFHandlerType::face_iterator cell_face = cell->face(face);
-              if (! cell->at_boundary(face) )
+              const bool periodic_neighbor = cell->has_periodic_neighbor(face);
+              if (! cell->at_boundary(face) || periodic_neighbor)
                 {
-                  typename DoFHandlerType::level_cell_iterator neighbor = cell->neighbor(face);
+                  typename DoFHandlerType::level_cell_iterator neighbor
+                    = cell->neighbor_or_periodic_neighbor(face);
 
                   // in 1d, we do not need to worry whether the neighbor
                   // might have children and then loop over those children.
-                  // rather, we may as well go straight to to cell behind
+                  // rather, we may as well go straight to the cell behind
                   // this particular cell's most terminal child
                   if (DoFHandlerType::dimension==1)
                     while (neighbor->has_children())
@@ -571,9 +573,10 @@ namespace DoFTools
                            sub_nr != cell_face->number_of_children();
                            ++sub_nr)
                         {
-                          const typename DoFHandlerType::level_cell_iterator
-                          sub_neighbor
-                            = cell->neighbor_child_on_subface (face, sub_nr);
+                          const typename DoFHandlerType::level_cell_iterator sub_neighbor
+                            = periodic_neighbor?
+                              cell->periodic_neighbor_child_on_subface (face, sub_nr):
+                              cell->neighbor_child_on_subface (face, sub_nr);
 
                           const unsigned int n_dofs_on_neighbor
                             = sub_neighbor->get_fe().dofs_per_cell;
@@ -598,9 +601,10 @@ namespace DoFTools
                     {
                       // Refinement edges are taken care of by coarser
                       // cells
-                      if (cell->neighbor_is_coarser(face) &&
-                          neighbor->subdomain_id() == cell->subdomain_id())
-                        continue;
+                      if ((!periodic_neighbor && cell->neighbor_is_coarser(face)) ||
+                          (periodic_neighbor && cell->periodic_neighbor_is_coarser(face)))
+                        if (neighbor->subdomain_id() == cell->subdomain_id())
+                          continue;
 
                       const unsigned int n_dofs_on_neighbor
                         = neighbor->get_fe().dofs_per_cell;
@@ -616,7 +620,7 @@ namespace DoFTools
                       // is not locally owned - otherwise, we touch each
                       // face twice and hence put the indices the other way
                       // around
-                      if (!cell->neighbor(face)->active()
+                      if (!cell->neighbor_or_periodic_neighbor(face)->active()
                           ||
                           (neighbor->subdomain_id() != cell->subdomain_id()))
                         {
@@ -755,7 +759,9 @@ namespace DoFTools
                   const typename DoFHandlerType::face_iterator
                   cell_face = cell->face (face_n);
 
-                  if (cell->at_boundary (face_n))
+                  const bool periodic_neighbor = cell->has_periodic_neighbor (face_n);
+
+                  if (cell->at_boundary (face_n) && (!periodic_neighbor))
                     {
                       for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
                         {
@@ -779,7 +785,7 @@ namespace DoFTools
                   else
                     {
                       typename DoFHandlerType::level_cell_iterator
-                      neighbor = cell->neighbor (face_n);
+                      neighbor = cell->neighbor_or_periodic_neighbor (face_n);
                       // If the cells are on the same level then only add to
                       // the sparsity pattern if the current cell is 'greater'
                       // in the total ordering.
@@ -796,7 +802,8 @@ namespace DoFTools
                       // if the neighbor is coarser than the current cell.
                       if (neighbor->level () != cell->level ()
                           &&
-                          !cell->neighbor_is_coarser (face_n))
+                          ((!periodic_neighbor && !cell->neighbor_is_coarser (face_n))
+                           ||(periodic_neighbor && !cell->periodic_neighbor_is_coarser (face_n))))
                         continue; // (the neighbor is finer)
 
                       const unsigned int
@@ -810,7 +817,9 @@ namespace DoFTools
                             {
                               const typename DoFHandlerType::level_cell_iterator
                               sub_neighbor
-                                = cell->neighbor_child_on_subface (face_n, sub_nr);
+                                = periodic_neighbor?
+                                  cell->periodic_neighbor_child_on_subface (face_n, sub_nr):
+                                  cell->neighbor_child_on_subface (face_n, sub_nr);
 
                               sub_neighbor->get_dof_indices (dofs_on_other_cell);
                               for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
@@ -1002,7 +1011,9 @@ namespace DoFTools
                 if (cell_face->user_flag_set ())
                   continue;
 
-                if (cell->at_boundary (face) )
+                const bool periodic_neighbor = cell->has_periodic_neighbor (face);
+
+                if (cell->at_boundary (face) && (!periodic_neighbor))
                   {
                     for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
                       for (unsigned int j=0; j<cell->get_fe().dofs_per_cell; ++j)
@@ -1019,14 +1030,17 @@ namespace DoFTools
                 else
                   {
                     typename dealii::hp::DoFHandler<dim,spacedim>::level_cell_iterator
-                    neighbor = cell->neighbor(face);
+                    neighbor = cell->neighbor_or_periodic_neighbor(face);
 
                     // Refinement edges are taken care of by coarser cells
-                    if (cell->neighbor_is_coarser(face))
+                    if ((!periodic_neighbor && cell->neighbor_is_coarser(face)) ||
+                        (periodic_neighbor && cell->periodic_neighbor_is_coarser(face)))
                       continue;
 
                     const unsigned int
-                    neighbor_face = cell->neighbor_of_neighbor(face);
+                    neighbor_face = periodic_neighbor?
+                                    cell->periodic_neighbor_of_periodic_neighbor(face):
+                                    cell->neighbor_of_neighbor(face);
 
                     if (cell_face->has_children())
                       {
@@ -1036,7 +1050,9 @@ namespace DoFTools
                           {
                             const typename dealii::hp::DoFHandler<dim,spacedim>::level_cell_iterator
                             sub_neighbor
-                              = cell->neighbor_child_on_subface (face, sub_nr);
+                              = periodic_neighbor?
+                                cell->periodic_neighbor_child_on_subface (face, sub_nr):
+                                cell->neighbor_child_on_subface (face, sub_nr);
 
                             dofs_on_other_cell.resize (sub_neighbor->get_fe().dofs_per_cell);
                             sub_neighbor->get_dof_indices (dofs_on_other_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.