]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
make patches in DoFTools
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 28 Sep 2011 05:14:05 +0000 (05:14 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 28 Sep 2011 05:14:05 +0000 (05:14 +0000)
git-svn-id: https://svn.dealii.org/trunk@24456 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 ddb98626acdb09bcf64974e8cdfe9c4b7d3935af..cd27f3b21d08ab52de1713ef9afbb6591e5bf1d0 100644 (file)
@@ -19,6 +19,7 @@
 #include <deal.II/base/index_set.h>
 #include <deal.II/base/point.h>
 #include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/lac/sparsity_pattern.h>
 #include <deal.II/dofs/function_map.h>
 #include <deal.II/dofs/dof_handler.h>
 #include <deal.II/fe/fe.h>
@@ -221,11 +222,19 @@ namespace DoFTools
                                   /**
                                    * Maximal number of degrees of
                                    * freedom on a cell.
+                                   *
+                                   * @relates DoFHandler
                                    */
   template <int dim, int spacedim>
   unsigned int
   max_dofs_per_cell (const DoFHandler<dim,spacedim> &dh);
 
+                                  /**
+                                   * Maximal number of degrees of
+                                   * freedom on a cell.
+                                   *
+                                   * @relates hp::DoFHandler
+                                   */
   template <int dim, int spacedim>
   unsigned int
   max_dofs_per_cell (const hp::DoFHandler<dim,spacedim> &dh);
@@ -239,6 +248,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates DoFHandler
                                    */
   template <int dim, int spacedim>
   unsigned int
@@ -252,6 +263,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates hp::DoFHandler
                                    */
   template <int dim, int spacedim>
   unsigned int
@@ -265,6 +278,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates DoFHandler
                                    */
   template <int dim, int spacedim>
   unsigned int
@@ -278,6 +293,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates hp::DoFHandler
                                    */
   template <int dim, int spacedim>
   unsigned int
@@ -292,6 +309,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates DoFHandler
                                    */
   template <int dim, int spacedim>
   unsigned int
@@ -306,6 +325,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates hp::DoFHandler
                                    */
   template <int dim, int spacedim>
   unsigned int
@@ -320,6 +341,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates DoFHandler
                                    */
   template <int dim, int spacedim>
   bool
@@ -334,6 +357,8 @@ namespace DoFTools
                                    * and hp DoFHandlers, to allow for a
                                    * uniform interface to query this
                                    * property.
+                                   *
+                                   * @relates hp::DoFHandler
                                    */
   template <int dim, int spacedim>
   bool
@@ -1455,6 +1480,21 @@ namespace DoFTools
                         const DH& dof_handler,
                         const unsigned int level);
   
+                                   /**
+                                    * Create a sparsity pattern which
+                                    * in each row lists the degrees of
+                                    * freedom associated to the
+                                    * corresponding cells
+                                    * corresponding to a vertex.
+                                    */
+   template <class DH>
+   void make_vertex_patches(SparsityPattern& block_list,
+                           const DH& dof_handler,
+                           const unsigned int level,
+                           const bool interior_dofs_only,
+                           const bool boundary_patches = true,
+                           const bool single_cell_patches = false);
+  
                                   //@}
                                   /**
                                    * Extract a vector that represents the
index de0e400b78ff6a00ce4f3027f940fa31fcf4a562..2d4eb64933e2aeb5532303dad292385c523bb0d5 100644 (file)
@@ -5859,9 +5859,146 @@ namespace DoFTools
            }
        }
   }
+
+  template <class DH, class Sparsity>
+  void make_cell_patches(Sparsity& block_list, const DH& dof_handler, const unsigned int level)
+  {
+    typename DH::cell_iterator cell;
+    typename DH::cell_iterator endc = dof_handler.end(level);
+    std::vector<unsigned int> indices;
+    
+    unsigned int i=0;
+    for (cell=dof_handler.begin(level); cell != endc; ++i, ++cell)
+      {
+       indices.resize(cell->get_fe().dofs_per_cell);
+       cell->get_mg_dof_indices(indices);
+       for (unsigned int j=0;j<indices.size();++j)
+         block_list.add(i,indices[j]);
+      }
+  }
+  
+  template <class DH>
+  void make_vertex_patches(
+    SparsityPattern& block_list,
+    const DH& dof_handler,
+    const unsigned int level,
+    const bool interior_only,
+    const bool boundary_patches,
+    const bool single_cell_patches)
+  {
+    typename DH::cell_iterator cell;
+    typename DH::cell_iterator endc = dof_handler.end(level);
+
+                                    // Vector mapping from vertex
+                                    // index in the triangulation to
+                                    // consecutive block indices on
+                                    // this level
+                                    // The number of cells at a vertex
+    std::vector<unsigned int> vertex_cell_count(dof_handler.get_tria().n_vertices(), 0);
+                                    // Is a vertex at the boundary?
+    std::vector<bool> vertex_boundary(dof_handler.get_tria().n_vertices(), false);
+    std::vector<unsigned int> vertex_mapping(dof_handler.get_tria().n_vertices(),
+                                            numbers::invalid_unsigned_int);
+                                    // Estimate for the number of
+                                    // dofs at this point
+    std::vector<unsigned int> vertex_dof_count(dof_handler.get_tria().n_vertices(), 0);
+    
+                                    // Identify all vertices active
+                                    // on this level and remember
+                                    // some data about them
+    for (cell=dof_handler.begin(level); cell != endc; ++cell)
+      for (unsigned int v=0;v<GeometryInfo<DH::dimension>::vertices_per_cell;++v)
+       {
+         const unsigned int vg = cell->vertex_index(v);
+         vertex_dof_count[vg] += cell->get_fe().dofs_per_cell;
+         ++vertex_cell_count[vg];
+         for (unsigned int d=0;d<DH::dimension;++d)
+           if (cell->at_boundary(GeometryInfo<DH::dimension>::vertex_to_face[v][d]))
+             vertex_boundary[vg] = true;
+       }
+                                    // From now on, onlly vertices
+                                    // with positive dof count are
+                                    // "in".
+    
+                                    // Remove vertices at boundaries
+                                    // or in corners
+    for (unsigned int vg=0;vg<vertex_dof_count.size();++vg)
+      if ((!single_cell_patches && vertex_cell_count[vg] < 2)
+         ||
+         (!boundary_patches && vertex_boundary[vg]))
+       vertex_dof_count[vg] = 0;
+
+                                    // Create a mapping from all
+                                    // vertices to the ones used here
+    unsigned int i=0;
+    for (unsigned int vg=0;vg<vertex_mapping.size();++vg)
+      if (vertex_dof_count[vg] != 0)
+       vertex_mapping[vg] = i++;
+
+                                    // Compactify dof count
+    for (unsigned int vg=0;vg<vertex_mapping.size();++vg)
+      if (vertex_dof_count[vg] != 0)
+       vertex_dof_count[vertex_mapping[vg]] = vertex_dof_count[vg];
+    
+                                    // Now that we have all the data,
+                                    // we reduce it to the part we
+                                    // actually want
+    vertex_dof_count.resize(i);
+    
+                                    // At this point, the list of
+                                    // patches is ready. Now we enter
+                                    // the dofs into the sparsity
+                                    // pattern.
+    block_list.reinit(vertex_dof_count.size(), dof_handler.n_dofs(level), vertex_dof_count);
+    
+    std::vector<unsigned int> indices;
+    std::vector<bool> exclude;
+    
+    for (cell=dof_handler.begin(level); cell != endc; ++cell)
+      {
+       const FiniteElement<DH::dimension>& fe = cell->get_fe();
+       indices.resize(fe.dofs_per_cell);
+       cell->get_mg_dof_indices(indices);
+       
+       for (unsigned int v=0;v<GeometryInfo<DH::dimension>::vertices_per_cell;++v)
+         {
+           const unsigned int vg = cell->vertex_index(v);
+           const unsigned int block = vertex_mapping[vg];
+           if (block == numbers::invalid_unsigned_int)
+             continue;
+           
+           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 d=0;d<DH::dimension;++d)
+                 {
+                   const unsigned int a_face = GeometryInfo<DH::dimension>::vertex_to_face[v][d];
+                   const unsigned int face = GeometryInfo<DH::dimension>::opposite_face[a_face];
+                   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(block, indices[j]);
+             }
+           else
+             {
+               for (unsigned int j=0;j<indices.size();++j)
+                 block_list.add(block, indices[j]);
+             }
+         }
+      }
+  }
 }
 
 
+
 // explicit instantiations
 
 #include "dof_tools.inst"
index 28287e1cf924f88742c2ab84003b38a6bd4a993b..ab1df6b2e5b1f77abbc12954b089b43d3e2da6b9 100644 (file)
@@ -113,7 +113,11 @@ for (SP : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS)
     DoFTools::make_flux_sparsity_pattern<hp::DoFHandler<deal_II_dimension>,SP>
     (const hp::DoFHandler<deal_II_dimension> &dof,
      SP    &sparsity);
-
+     
+    template void
+      DoFTools::make_cell_patches<MGDoFHandler<deal_II_dimension>,SP>
+      (SP&, const MGDoFHandler<deal_II_dimension>&, unsigned int);
+    
 #if deal_II_dimension > 1
 
     template void
@@ -236,6 +240,9 @@ for (DH : DOFHANDLERS; deal_II_dimension : DIMENSIONS)
 
 for (deal_II_dimension : DIMENSIONS)
 {
+  template
+  void DoFTools::make_vertex_patches (SparsityPattern&, const MGDoFHandler<deal_II_dimension>&,
+  unsigned int, bool, bool, bool);
 
 #if deal_II_dimension < 3
 template

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.