From: leicht Date: Wed, 13 Dec 2006 12:48:10 +0000 (+0000) Subject: new function: delete_duplicated_vertices X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f70db9bb22cc77d2de9151e71d38b8fbf4d79ebc;p=dealii-svn.git new function: delete_duplicated_vertices git-svn-id: https://svn.dealii.org/trunk@14233 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/grid/grid_tools.h b/deal.II/deal.II/include/grid/grid_tools.h index 60dbd64b42..d8e10df557 100644 --- a/deal.II/deal.II/include/grid/grid_tools.h +++ b/deal.II/deal.II/include/grid/grid_tools.h @@ -115,6 +115,37 @@ class GridTools std::vector > &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 + * GridIn::read_* 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 + static + void delete_duplicated_vertices (std::vector > &all_vertices, + std::vector > &cells, + SubCellData &subcelldata, + std::vector &considered_vertices, + double tol=1e-12); + /** * Transform the vertices of the * given triangulation by diff --git a/deal.II/deal.II/source/grid/grid_tools.cc b/deal.II/deal.II/source/grid/grid_tools.cc index 1d1b50d3ff..fb0171b79e 100644 --- a/deal.II/deal.II/source/grid/grid_tools.cc +++ b/deal.II/deal.II/source/grid/grid_tools.cc @@ -245,7 +245,61 @@ GridTools::cell_measure (const std::vector > &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 > &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 double @@ -256,7 +310,6 @@ GridTools::cell_measure(const std::vector > &all_vertices, return 0.; } -#endif template @@ -315,6 +368,75 @@ GridTools::delete_unused_vertices (std::vector > &vertices, +template +void +GridTools::delete_duplicated_vertices (std::vector > &vertices, + std::vector > &cells, + SubCellData &subcelldata, + std::vector &considered_vertices, + double tol) +{ + // create a vector of vertex + // indices. initialize it to the identity, + // later on change that if possible. + std::vector new_vertex_numbers(vertices.size()); + //TODO: find a better way to initialize this + //vector + for (unsigned int i=0; i::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 > & std::vector > &, SubCellData &); +template +void GridTools::delete_duplicated_vertices (std::vector > &, + std::vector > &, + SubCellData &, + std::vector &, + double); + template void GridTools::shift (const Point &, Triangulation &); diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index e0b9c9c9d2..e564a39fd6 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -397,6 +397,12 @@ inconvenience this causes.

base

    +
  1. Fixed: If no substring is found with a width smaller than the given + threshold, the Utilities::break_text_into_lines function now + returns the smallest possible substring (larger than the threshold). +
    (Tobias Leicht 2006/12/13) +

    +
  2. Changed: We used to set preprocessor variables PACKAGE_NAME, PACKAGE_TARNAME, PACKAGE_VERSION, PACKAGE_STRING, @@ -922,6 +928,16 @@ inconvenience this causes.

      +
    1. New: There is now a function + GridTools::delete_duplicated_vertices 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. +
      + (Tobias Leicht 2006/12/13) +

      +
    2. New: There is now a function GridGenerator::subdivided_hyper_rectangle that produces a non-uniformly subdivided rectangle, ideally