]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement a second function for partitioning a triangulation.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 6 Aug 2008 17:15:08 +0000 (17:15 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 6 Aug 2008 17:15:08 +0000 (17:15 +0000)
git-svn-id: https://svn.dealii.org/trunk@16502 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_tools.h
deal.II/deal.II/source/grid/grid_tools.cc
deal.II/doc/news/changes.h

index a4bce05c5d906049134ea3f8b7458a3a4b4d51ea..9e6a88e70e9b74bcf6ea1935b9d39410b7bb971e 100644 (file)
@@ -33,6 +33,8 @@ namespace hp
   template <int> class MappingCollection;
 }
 
+class SparsityPattern;
+
 
 /**
  * This class is a collection of algorithms working on triangulations,
@@ -446,7 +448,7 @@ class GridTools
                                       * values between zero and
                                       * @p n_partitions-1. You can access the
                                       * subdomain id of a cell by using
-                                      * <tt>cell->subdomain_id()</tt>.
+                                      * <tt>cell-@>subdomain_id()</tt>.
                                       *
                                       * This function will generate an error
                                       * if METIS is not installed unless
@@ -461,6 +463,88 @@ class GridTools
     static
     void partition_triangulation (const unsigned int  n_partitions,
                                   Triangulation<dim> &triangulation);
+
+                                    /**
+                                     * This function does the same as the
+                                     * previous one, i.e. it partitions a
+                                     * triangulation using METIS into a
+                                     * number of subdomains identified by the
+                                     * <code>cell-@>subdomain_id()</code>
+                                     * flag.
+                                     *
+                                     * The difference to the previous
+                                     * function is the second argument, a
+                                     * sparsity pattern that represents the
+                                     * connectivity pattern between cells.
+                                     *
+                                     * While the function above builds it
+                                     * directly from the triangulation by
+                                     * considering which cells neighbor each
+                                     * other, this function can take a more
+                                     * refined connectivity graph. The
+                                     * sparsity pattern needs to be of size
+                                     * $N\times N$, where $N$ is the number
+                                     * of active cells in the
+                                     * triangulation. If the sparsity pattern
+                                     * contains an entry at position $(i,j)$,
+                                     * then this means that cells $i$ and $j$
+                                     * (in the order in which they are
+                                     * traversed by active cell iterators)
+                                     * are to be considered connected; METIS
+                                     * will then try to partition the domain
+                                     * in such a way that (i) the subdomains
+                                     * are of roughly equal size, and (ii) a
+                                     * minimal number of connections are
+                                     * broken.
+                                     *
+                                     * This function is mainly useful in
+                                     * cases where connections between cells
+                                     * exist that are not present in the
+                                     * triangulation alone (otherwise the
+                                     * previous function would be the simpler
+                                     * one to use). Such connections may
+                                     * include that certain parts of the
+                                     * boundary of a domain are coupled
+                                     * through symmetric boundary conditions
+                                     * or integrals (e.g. friction contact
+                                     * between the two sides of a crack in
+                                     * the domain), or if a numerical scheme
+                                     * is used that not only connects
+                                     * immediate neighbors but a larger
+                                     * neighborhood of cells (e.g. when
+                                     * solving integral equations).
+                                     *
+                                     * In addition, this function may be
+                                     * useful in cases where the default
+                                     * sparsity pattern is not entirely
+                                     * sufficient. This can happen because
+                                     * the default is to just consider face
+                                     * neighbors, not neighboring cells that
+                                     * are connected by edges or
+                                     * vertices. While the latter couple when
+                                     * using continuous finite elements, they
+                                     * are typically still closely connected
+                                     * in the neighborship graph, and METIS
+                                     * will not usually cut important
+                                     * connections in this case. However, if
+                                     * there are vertices in the mesh where
+                                     * many cells (many more than the common
+                                     * 4 or 6 in 2d and 3d, respectively)
+                                     * come together, then there will be a
+                                     * significant number of cells that are
+                                     * connected across a vertex, but several
+                                     * degrees removed in the connectivity
+                                     * graph built only using face
+                                     * neighbors. In a case like this, METIS
+                                     * may sometimes make bad decisions and
+                                     * you may want to build your own
+                                     * connectivity graph.
+                                     */
+    template <int dim>
+    static
+    void partition_triangulation (const unsigned int     n_partitions,
+                                 const SparsityPattern &cell_connection_graph,
+                                  Triangulation<dim>    &triangulation);
     
                                      /**
                                       * For each active cell, return in the
index c552eaf894116bf624b172ab18ecccbc89ca602c..868ea6c96458164e79027d76f30b607dfca82fb2 100644 (file)
@@ -839,61 +839,94 @@ partition_triangulation (const unsigned int  n_partitions,
                                    // we decompose the domain by first
                                    // generating the connection graph of all
                                    // cells with their neighbors, and then
-                                   // passing this graph off to METIS. To make
-                                   // things a little simpler and more
-                                   // general, we let the function
-                                   // DoFTools:make_flux_sparsity_pattern
-                                   // function generate the connection graph
-                                   // for us and reuse the SparsityPattern
-                                   // data structure for the connection
-                                   // graph. The connection structure of the
-                                   // mesh is obtained by using a fake
-                                   // piecewise constant finite element
-                                   //
-                                   // Since in 3d the generation of a
-                                   // sparsity pattern can be expensive, we
-                                   // take the detour of the compressed
-                                   // sparsity pattern, which is a little
-                                   // slower but more efficient in terms of
-                                   // memory
-  FE_DGQ<dim> fake_q0(0);
-  DoFHandler<dim> dof_handler (triangulation);
-  dof_handler.distribute_dofs (fake_q0);
-  Assert (dof_handler.n_dofs() == triangulation.n_active_cells(),
-          ExcInternalError());
+                                   // passing this graph off to METIS.
+                                  //
+                                  // as built in this function, we only
+                                  // consider face neighbors, which leads to
+                                  // a fixed number of entries per row (don't
+                                  // forget that each cell couples with
+                                  // itself)
+  SparsityPattern cell_connectivity (triangulation.n_active_cells(),
+                                    triangulation.n_active_cells(),
+                                    GeometryInfo<dim>::faces_per_cell+1);
+
+                                  // next we have to build a mapping from the
+                                  // list of cells to their indices. for
+                                  // this, use the user_index field
+  std::vector<unsigned int> saved_user_indices;
+  triangulation.save_user_indices (saved_user_indices);
+  unsigned int index = 0;
+  for (typename Triangulation<dim>::active_cell_iterator
+         cell = triangulation.begin_active();
+       cell != triangulation.end(); ++cell, ++index)
+    cell->set_user_index (index);
+
+                                  // next loop over all cells and their
+                                  // neighbors to build the sparsity pattern
+  index = 0;
+  for (typename Triangulation<dim>::active_cell_iterator
+         cell = triangulation.begin_active();
+       cell != triangulation.end(); ++cell, ++index)
+    {
+      cell_connectivity.add (index, index);
+      for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+       if (cell->at_boundary(f) == false)
+         cell_connectivity.add (index,
+                                cell->neighbor(f)->user_index());
+    }
   
-  CompressedSparsityPattern csp (dof_handler.n_dofs(),
-                                 dof_handler.n_dofs());
-  DoFTools::make_flux_sparsity_pattern (dof_handler, csp);
+                                  // now compress the so-built connectivity
+                                  // pattern and restore user indices
+  cell_connectivity.compress ();
+  triangulation.save_user_indices (saved_user_indices);
   
-  SparsityPattern sparsity_pattern;
-  sparsity_pattern.copy_from (csp);
+                                  // finally defer to the other function for
+                                  // partitioning and assigning subdomain ids
+  partition_triangulation (n_partitions,
+                          cell_connectivity,
+                          triangulation);
+}
+
+
+
+template <int dim>
+void
+GridTools::
+partition_triangulation (const unsigned int     n_partitions,
+                        const SparsityPattern &cell_connection_graph,
+                         Triangulation<dim>    &triangulation)
+{
+  Assert (n_partitions > 0, ExcInvalidNumberOfPartitions(n_partitions));
+  Assert (cell_connection_graph.n_rows() == triangulation.n_active_cells(),
+         ExcMessage ("Connectivity graph has wrong size"));
+  Assert (cell_connection_graph.n_cols() == triangulation.n_active_cells(),
+         ExcMessage ("Connectivity graph has wrong size"));  
+  
+                                   // check for an easy return
+  if (n_partitions == 1)
+    {
+      for (typename Triangulation<dim>::active_cell_iterator
+             cell = triangulation.begin_active();
+           cell != triangulation.end(); ++cell)
+        cell->set_subdomain_id (0);
+      return;
+    }
 
                                    // partition this connection graph and get
                                    // back a vector of indices, one per degree
                                    // of freedom (which is associated with a
                                    // cell)
   std::vector<unsigned int> partition_indices (triangulation.n_active_cells());
-  sparsity_pattern.partition (n_partitions,  partition_indices);
+  cell_connection_graph.partition (n_partitions,  partition_indices);
 
                                    // finally loop over all cells and set the
-                                   // subdomain ids. for this, get the DoF
-                                   // index of each cell and extract the
-                                   // subdomain id from the vector obtained
-                                   // above
+                                   // subdomain ids
   std::vector<unsigned int> dof_indices(1);
-  for (typename DoFHandler<dim>::active_cell_iterator
-         cell = dof_handler.begin_active();
-       cell != dof_handler.end(); ++cell)
-    {
-      cell->get_dof_indices(dof_indices);
-      Assert (dof_indices[0] < triangulation.n_active_cells(),
-              ExcInternalError());
-      Assert (partition_indices[dof_indices[0]] < n_partitions,
-              ExcInternalError());
-      
-      cell->set_subdomain_id (partition_indices[dof_indices[0]]);
-    }
+  unsigned int index = 0;
+  for (typename Triangulation<dim>::active_cell_iterator
+         cell = triangulation.begin_active();
+       cell != triangulation.end(); ++cell, ++index)
+    cell->set_subdomain_id (partition_indices[index]);
 }
 
 
@@ -1278,6 +1311,12 @@ void
 GridTools::partition_triangulation (const unsigned int,
                                     Triangulation<deal_II_dimension> &);
 
+template
+void
+GridTools::partition_triangulation (const unsigned int,
+                                   const SparsityPattern &,
+                                    Triangulation<deal_II_dimension> &);
+
 template
 void
 GridTools::
index 0394b0109d552302052d1207a495e69dc636350f..eda6995d4f3f4ae458d97d42cb049bb5a5e5d3a3 100644 (file)
@@ -234,18 +234,28 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li>
+  <p>
+  New: There is a second GridTools::partition_triangulation
+  function that takes a cell connectivity pattern as argument, rather
+  than computing it itself as the existing function. Use cases are
+  discussed in the documentation of the new function.
+  (WB 2008/08/06)
+  <li>
+
   <li>
   <p>
   Fixed: GridTools::find_cells_adjacent_to_vertex had a bug that 
   prevented its correct functioning in three dimensions. Some 
   cases were left out due to uncorrect assumptions on the various 
   refinement possibilities.
+  <br>
   (Luca Heltai 2008/07/17)
   <li>
 
   <p>
   New: There is now a new
-  Triangulation::prevent_distorted_boundary_cells function which
+  Triangulation::prevent_distorted_boundary_cells function which is
   only useful in case of anisotropic refinement. At the boundary
   of the domain, the new point on the face may be far inside the
   current cell, if the boundary has a strong curvature. If we

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.