]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
new patches
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 30 Oct 2011 19:06:25 +0000 (19:06 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 30 Oct 2011 19:06:25 +0000 (19:06 +0000)
git-svn-id: https://svn.dealii.org/trunk@24703 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/dofs/dof_tools.h
deal.II/source/dofs/dof_tools.cc
deal.II/source/dofs/dof_tools.inst.in

index e387e694c48b8d6a8a063adfc714511f1d3ac2ac..540158a2d7a616e6e8a2d955760ee2e218584c94 100644 (file)
@@ -1492,10 +1492,14 @@ namespace DoFTools
   
                                    /**
                                     * Create a sparsity pattern which
-                                    * in each row lists the degrees of
-                                    * freedom associated to the
-                                    * corresponding cells
-                                    * corresponding to a vertex.
+                                    * in each row lists the degrees
+                                    * of freedom associated to the
+                                    * cells corresponding to a
+                                    * vertex. The sparsity pattern
+                                    * may be empty when entering this
+                                    * function and will be
+                                    * reinitialized to the correct
+                                    * size.
                                     *
                                     * The function has some boolean
                                     * arguments (lsited below)
@@ -1553,6 +1557,99 @@ namespace DoFTools
                            const bool level_boundary_patches = false,
                            const bool single_cell_patches = false);
   
+                                   /**
+                                    * Create a sparsity pattern which
+                                    * in each row lists the degrees of
+                                    * freedom associated to the
+                                    * cells which are the children of
+                                    * the same cell. The
+                                    * sparsity pattern may be empty
+                                    * when entering this function and
+                                    * will be reinitialized to the
+                                    * correct size.
+                                    *
+                                    * The function has some boolean
+                                    * arguments (lsited below)
+                                    * controlling details of the
+                                    * generated patches. The default
+                                    * settings are those for
+                                    * Arnold-Falk-Winther type
+                                    * smoothers for divergence and
+                                    * curl conforming finite elements
+                                    * with essential boundary
+                                    * conditions. Other applications
+                                    * are possible, in particular
+                                    * changing
+                                    * <tt>boundary_dofs</tt> for
+                                    * non-essential boundary
+                                    * conditions.
+                                    *
+                                    * Since the patches are defined
+                                    * through refinement, th 
+                                    *
+                                    * @arg <tt>block_list</tt>: the
+                                    * SparsityPattern into which the
+                                    * patches will be stored.
+                                    * @arg <tt>dof_handler</tt>: The
+                                    * multilevel dof handler
+                                    * providing the topology operated
+                                    * on.
+                                    * @arg
+                                    * <tt>interior_dofs_only</tt>:
+                                    * for each patch of cells around
+                                    * a vertex, collect only the
+                                    * interior degrees of freedom of
+                                    * the patch and disregard those
+                                    * on the boundary of the
+                                    * patch. This is for instance the
+                                    * setting for smoothers of
+                                    * Arnold-Falk-Winther type.
+                                    * @arg <tt>boundary_dofs</tt>:
+                                    * include degrees of freedom,
+                                    * which would have excluded by
+                                    * <tt>interior_dofs_only</tt>,
+                                    * but are lying on the boundary
+                                    * of the domain, and thus need smoothing.
+                                    */
+   template <class DH>
+   void make_child_patches(SparsityPattern& block_list,
+                           const DH& dof_handler,
+                           const unsigned int level,
+                           const bool interior_dofs_only,
+                           const bool boundary_dofs = false);
+
+                                  /**
+                                   * Create a block list with only a
+                                   * single patch, which in turn
+                                   * contains all degrees of freedom
+                                   * on the given level.
+                                   *
+                                   * This function is mostly a
+                                   * closure on level 0 for functions
+                                   * like make_child_patches() and
+                                   * make_vertex_patches(), which may
+                                   * produce an empty patch list.
+                                   *
+                                   * @arg <tt>block_list</tt>: the
+                                   * SparsityPattern into which the
+                                   * patches will be stored.
+                                   * @arg <tt>dof_handler</tt>: The
+                                   * multilevel dof handler
+                                   * providing the topology operated
+                                   * on.
+                                   * @arg <tt>level</tt> The grid
+                                   * level used for building the list.
+                                   * @arg
+                                   * <tt>interior_dofs_only</tt>:
+                                   * if true, exclude degrees of freedom on
+                                   * the boundary of the domain.
+                                   */
+  template <class DH>
+   void make_single_patch(SparsityPattern& block_list,
+                         const DH& dof_handler,
+                         const unsigned int level,
+                         const bool interior_dofs_only = false);
+  
                                   //@}
                                   /**
                                    * Extract a vector that represents the
index 2de1b6646ef735228a6e34252d27a4e86c66fee7..2a89f8b8c4fa71ab4d7932e8a0f95a0beaca626c 100644 (file)
@@ -5877,6 +5877,100 @@ namespace DoFTools
       }
   }
   
+  
+  template <class DH>
+  void make_single_patch(
+    SparsityPattern& block_list,
+    const DH& dof_handler,
+    const unsigned int level,
+    const bool interior_only)
+    {
+      const FiniteElement<DH::dimension>& fe = dof_handler.get_fe();
+      block_list.reinit(1, dof_handler.n_dofs(level), dof_handler.n_dofs(level));
+      typename DH::cell_iterator cell;
+      typename DH::cell_iterator endc = dof_handler.end(level);
+      
+      std::vector<unsigned int> indices;
+      std::vector<bool> exclude;
+      
+      for (cell=dof_handler.begin(level); cell != endc;++cell)
+      {
+       indices.resize(cell->get_fe().dofs_per_cell);
+       cell->get_mg_dof_indices(indices);
+       
+       if (interior_only)
+         {
+                                            // Exclude degrees of
+                                            // freedom on faces
+                                            // opposite to the vertex
+           exclude.resize(fe.dofs_per_cell);
+           std::fill(exclude.begin(), exclude.end(), false);
+           const unsigned int dpf = fe.dofs_per_face;
+           
+           for (unsigned int face=0;face<GeometryInfo<DH::dimension>::faces_per_cell;++face)
+             if (cell->at_boundary(face) || cell->neighbor(face)->level() != cell->level())
+               for (unsigned int i=0;i<dpf;++i)
+                 exclude[fe.face_to_cell_index(i,face)] = true;
+           for (unsigned int j=0;j<indices.size();++j)
+             if (!exclude[j])
+               block_list.add(0, indices[j]);
+         }
+       else
+         {
+           for (unsigned int j=0;j<indices.size();++j)
+             block_list.add(0, indices[j]);
+         }
+      }
+    }
+
+  
+  template <class DH>
+  void make_child_patches(
+    SparsityPattern& block_list,
+    const DH& dof_handler,
+    const unsigned int level,
+    const bool interior_dofs_only,
+    const bool boundary_dofs)
+  {
+    Assert(level > 0 && level < dof_handler.get_tria().n_levels(),
+          ExcIndexRange(level, 1, dof_handler.get_tria().n_levels()));
+
+    typename DH::cell_iterator cell;
+    typename DH::cell_iterator pcell = dof_handler.begin(level-1);
+    typename DH::cell_iterator endc = dof_handler.end(level-1);
+    
+    std::vector<unsigned int> indices;
+    std::vector<bool> exclude;
+    
+    for (unsigned int block = 0;pcell != endc;++pcell)
+      {
+       if (pcell->has_children()) continue;
+       for (unsigned int child=0;child<pcell->n_children();++child)
+         {
+           cell = pcell->child(child);
+                                            // For hp, only this line
+                                            // here would have to be
+                                            // replaced.
+           const FiniteElement<DH::dimension>& fe = dof_handler.get_fe();
+           const unsigned int n_dofs = fe.dofs_per_cell;
+           indices.resize(n_dofs);
+           exclude.resize(n_dofs);
+           std::fill(exclude.begin(), exclude.end(), false);
+           cell->get_mg_dof_indices(indices);
+
+           Assert(!interior_dofs_only, ExcNotImplemented());
+           Assert(!boundary_dofs, ExcNotImplemented());
+           
+           for (unsigned int i=0;i<n_dofs;++i)
+             if (!exclude[i])
+               block_list.add(block, indices[i]);
+         }
+       ++block;
+      }
+  }
+  
+  
+    
   template <class DH>
   void make_vertex_patches(
     SparsityPattern& block_list,
index d1d29ef355c9930b4d6eeda2245f8d26466dddbf..c70c279e45de0404083e4a80b12e67f695ecaa9f 100644 (file)
@@ -244,6 +244,14 @@ for (deal_II_dimension : DIMENSIONS)
   void DoFTools::make_vertex_patches (SparsityPattern&, const MGDoFHandler<deal_II_dimension>&,
   unsigned int, bool, bool, bool, bool);
 
+  template
+  void DoFTools::make_single_patch (SparsityPattern&, const MGDoFHandler<deal_II_dimension>&,
+  unsigned int, bool);
+  
+  template
+  void DoFTools::make_child_patches(SparsityPattern&, const MGDoFHandler<deal_II_dimension>&,
+  unsigned int, bool, bool);
+
 #if deal_II_dimension < 3
 template
 void

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.