]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
new function: delete_duplicated_vertices
authorleicht <leicht@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 13 Dec 2006 12:48:10 +0000 (12:48 +0000)
committerleicht <leicht@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 13 Dec 2006 12:48:10 +0000 (12:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@14233 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.html

index 60dbd64b4222d9698d4039095bc4cc3759030284..d8e10df5577d800b967e99700ee30f3e9d1538c0 100644 (file)
@@ -115,6 +115,37 @@ class GridTools
                                 std::vector<CellData<dim> > &cells,
                                 SubCellData                 &subcelldata);
     
+                                    /**
+                                     * Remove vertices that are duplicated,
+                                     * due to the input of a structured grid,
+                                     * for example. If these vertices are not
+                                     * removed, the faces bounded by these
+                                     * vertices become part of the boundary,
+                                     * even if they are in the interior of
+                                     * the mesh.
+                                     *
+                                     * This function is called by some
+                                     * <tt>GridIn::read_*</tt> functions. Only
+                                     * the vertices with indices in @p
+                                     * considered_vertices are tested for
+                                     * equality. This speeds up the algorithm,
+                                     * which is quadratic and thus quite slow
+                                     * to begin with. However, if you whish to
+                                     * consider all vertices, simply pass an
+                                     * empty vector.
+                                     *
+                                     * Two vertices are considered equal, if
+                                     * their difference in each coordinate
+                                     * direction is less than @p tol.
+                                     */
+    template <int dim>
+    static
+    void delete_duplicated_vertices (std::vector<Point<dim> >    &all_vertices,
+                                    std::vector<CellData<dim> > &cells,
+                                    SubCellData                 &subcelldata,
+                                    std::vector<unsigned int>   &considered_vertices,
+                                    double                       tol=1e-12);
+    
                                     /**
                                      * Transform the vertices of the
                                      * given triangulation by
index 1d1b50d3ff9c69b1b6fd77aef3817d9724428de0..fb0171b79e54e1717cc96120d6cfcfe15906754d 100644 (file)
@@ -245,7 +245,61 @@ GridTools::cell_measure (const std::vector<Point<3> > &all_vertices,
   return (t34+t64+t95+t125+t156+t181+t207+t228)/12.;
 }
 
-#else
+#endif
+
+
+#if deal_II_dimension == 2
+
+template <>
+double
+GridTools::cell_measure(const std::vector<Point<2> > &all_vertices,
+                       const unsigned int (&vertex_indices) [GeometryInfo<2>::vertices_per_cell])
+{
+/*
+  Get the computation of the measure by this little Maple script. We
+  use the blinear mapping of the unit quad to the real quad. However,
+  every transformation mapping the unit faces to straight lines should
+  do.
+
+  Remember that the area of the quad is given by
+  \int_K 1 dx dy  = \int_{\hat K} |det J| d(xi) d(eta)
+
+  # x and y are arrays holding the x- and y-values of the four vertices
+  # of this cell in real space. 
+  x := array(0..3);
+  y := array(0..3);
+  tphi[0] := (1-xi)*(1-eta):
+  tphi[1] :=     xi*(1-eta):
+  tphi[2] := (1-xi)*eta:
+  tphi[3] :=     xi*eta:
+  x_real := sum(x[s]*tphi[s], s=0..3):
+  y_real := sum(y[s]*tphi[s], s=0..3):
+  detJ := diff(x_real,xi)*diff(y_real,eta) - diff(x_real,eta)*diff(y_real,xi):
+
+  measure := simplify ( int ( int (detJ, xi=0..1), eta=0..1)):
+  readlib(C):
+
+  C(measure, optimized);
+
+  additional optimizaton: divide by 2 only one time
+*/
+
+  const double x[4] = { all_vertices[vertex_indices[0]](0),
+                       all_vertices[vertex_indices[1]](0),
+                       all_vertices[vertex_indices[2]](0),
+                       all_vertices[vertex_indices[3]](0)};
+
+  const double y[4] = { all_vertices[vertex_indices[0]](1),
+                       all_vertices[vertex_indices[1]](1),
+                       all_vertices[vertex_indices[2]](1),
+                       all_vertices[vertex_indices[3]](1)};
+
+  return (-x[1]*y[0]+x[1]*y[3]+y[0]*x[2]+x[0]*y[1]-x[0]*y[2]-y[1]*x[3]-x[2]*y[3]+x[3]*y[2])/2;
+
+}
+
+#endif
+
 
 template <int dim>
 double
@@ -256,7 +310,6 @@ GridTools::cell_measure(const std::vector<Point<dim> > &all_vertices,
   return 0.;
 }
 
-#endif
 
 
 template <int dim>
@@ -315,6 +368,75 @@ GridTools::delete_unused_vertices (std::vector<Point<dim> >    &vertices,
 
 
 
+template <int dim>
+void
+GridTools::delete_duplicated_vertices (std::vector<Point<dim> >    &vertices,
+                                      std::vector<CellData<dim> > &cells,
+                                      SubCellData                 &subcelldata,
+                                      std::vector<unsigned int>   &considered_vertices,
+                                      double                       tol)
+{
+                                  // create a vector of vertex
+                                  // indices. initialize it to the identity,
+                                  // later on change that if possible.
+  std::vector<unsigned int> new_vertex_numbers(vertices.size());
+                                  //TODO: find a better way to initialize this
+                                  //vector
+  for (unsigned int i=0; i<vertices.size(); ++i)
+    new_vertex_numbers[i]=i;
+
+                                  // if the considered_vertices vector is
+                                  // empty, consider all vertices
+  if (considered_vertices.size()==0)
+    considered_vertices=new_vertex_numbers;
+
+                                  // now loop over all vertices to be
+                                  // considered and try to find an identical
+                                  // one
+  for (unsigned int i=0; i<considered_vertices.size(); ++i)
+    {
+      if (new_vertex_numbers[considered_vertices[i]]!=considered_vertices[i])
+                                        // this vertex has been identified with
+                                        // another one already, skip it in the
+                                        // test
+       continue;
+                                      // this vertex is not identified with
+                                      // another one so far. search in the list
+                                      // of remaining vertices. if a duplicate
+                                      // vertex is found, set the new vertex
+                                      // index for that vertex to this vertex'
+                                      // index.
+      for (unsigned int j=i+1; j<considered_vertices.size(); ++j)
+       {
+         bool equal=true;
+         for (unsigned int d=0; d<dim; ++d)
+           equal &= (fabs(vertices[considered_vertices[j]](d)-vertices[considered_vertices[i]](d))<tol);
+         if (equal)
+           {
+             new_vertex_numbers[considered_vertices[j]]=considered_vertices[i];
+             // we do not suppose, that there might be another duplicate
+             // vertex, so break here
+             break;
+           }
+       }               
+    }
+
+                                  // now we got a renumbering list. simply
+                                  // renumber all vertices (non-duplicate
+                                  // vertices get renumbered to themselves, so
+                                  // nothing bad happens). after that, the
+                                  // duplicate vertices will be unused, so call
+                                  // delete_unused_vertices() to do that part
+                                  // of the job.
+  for (unsigned int c=0; c<cells.size(); ++c)
+    for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+      cells[c].vertices[v]=new_vertex_numbers[cells[c].vertices[v]];
+
+  delete_unused_vertices(vertices,cells,subcelldata);
+}
+
+
+
 // define some transformations in an anonymous namespace
 #ifdef DEAL_II_ANON_NAMESPACE_BOGUS_WARNING
   namespace TRANS
@@ -967,6 +1089,13 @@ void GridTools::delete_unused_vertices (std::vector<Point<deal_II_dimension> > &
                                        std::vector<CellData<deal_II_dimension> > &,
                                        SubCellData &);
 
+template
+void GridTools::delete_duplicated_vertices (std::vector<Point<deal_II_dimension> > &,
+                                           std::vector<CellData<deal_II_dimension> > &,
+                                           SubCellData &,
+                                           std::vector<unsigned int> &,
+                                           double);
+
 template
 void GridTools::shift<deal_II_dimension> (const Point<deal_II_dimension> &,
                                          Triangulation<deal_II_dimension> &);
index e0b9c9c9d22a2bd07e159a017a732bc4e1fdc03a..e564a39fd61d9142ccd4e253fb5fcf4bc151a7c2 100644 (file)
@@ -397,6 +397,12 @@ inconvenience this causes.
 <h3>base</h3>
 
 <ol>
+  <li> <p> Fixed: If no substring is found with a width smaller than the given
+       threshold, the <code>Utilities::break_text_into_lines</code> function now
+       returns the smallest possible substring (larger than the threshold).
+       <br> (Tobias Leicht 2006/12/13)
+       </p>
+
   <li> <p> Changed: We used to set preprocessor variables
        <code>PACKAGE_NAME</code>, <code>PACKAGE_TARNAME</code>,
        <code>PACKAGE_VERSION</code>, <code>PACKAGE_STRING</code>,
@@ -922,6 +928,16 @@ inconvenience this causes.
 
 <ol>
 
+  <li> <p> New: There is now a function
+       <code>GridTools::delete_duplicated_vertices</code> that deletes duplicate
+       vertices which can occur if structured grids are read into deal, leading
+       to unwanted interior boundaries. In order to reduce the effort of the
+       quadratic algorithm, a list of vertices to be considered can be supplied
+       if available.  
+                <br>
+                (Tobias Leicht 2006/12/13)
+       </p>
+
   <li> <p> New: There is now a function <code
        class="member">GridGenerator::subdivided_hyper_rectangle</code>
        that produces a non-uniformly subdivided rectangle, ideally

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.