static void skip_comment_lines (std::istream &in,
const char comment_start);
- /**
- * Remove vertices that are not
- * referenced by any of the
- * cells. This function is called
- * by all <tt>read_*</tt> functions to
- * eliminate vertices that are
- * listed in the input files but
- * are not used by the cells in
- * the input file. While these
- * vertices should not be in the
- * input from the beginning, they
- * sometimes are, most often when
- * some cells have been removed
- * by hand without wanting to
- * update the vertex lists, as
- * they might be lengthy.
- *
- * This function is called by all
- * <tt>read_*</tt> functions as the
- * triangulation class requires
- * them to be called with used
- * vertices only. This is so,
- * since the vertices are copied
- * verbatim by that class, so we
- * have to eliminate unused
- * vertices beforehand.
- */
- static void delete_unused_vertices (std::vector<Point<dim> > &vertices,
- std::vector<CellData<dim> > &cells,
- SubCellData &subcelldata);
-
/**
* This function can write the
* raw cell data objects created
* individual functions for more information.
*
* @ingroup grid
- * @author Wolfgang Bangerth, 2001, 2003, 2004
+ * @author Wolfgang Bangerth, 2001, 2003, 2004, Ralf Hartmann, 2005
*/
class GridTools
{
double cell_measure(const std::vector<Point<dim> > &all_vertices,
const int vertex_indices[GeometryInfo<dim>::vertices_per_cell]);
+ /**
+ * Remove vertices that are not
+ * referenced by any of the
+ * cells. This function is called
+ * by all <tt>GridIn::read_*</tt>
+ * functions to eliminate
+ * vertices that are listed in
+ * the input files but are not
+ * used by the cells in the input
+ * file. While these vertices
+ * should not be in the input
+ * from the beginning, they
+ * sometimes are, most often when
+ * some cells have been removed
+ * by hand without wanting to
+ * update the vertex lists, as
+ * they might be lengthy.
+ *
+ * This function is called by all
+ * <tt>GridIn::read_*</tt>
+ * functions as the triangulation
+ * class requires them to be
+ * called with used vertices
+ * only. This is so, since the
+ * vertices are copied verbatim
+ * by that class, so we have to
+ * eliminate unused vertices
+ * beforehand.
+ */
+ template <int dim>
+ static
+ void delete_unused_vertices (std::vector<Point<dim> > &vertices,
+ std::vector<CellData<dim> > &cells,
+ SubCellData &subcelldata);
+
/**
* Transform the vertices of the
* given triangulation by
#include <grid/grid_in.h>
#include <grid/tria.h>
#include <grid/grid_reordering.h>
+#include <grid/grid_tools.h>
#include <map>
#include <algorithm>
AssertThrow (in, ExcIO());
// do some clean-up on vertices...
- delete_unused_vertices (vertices, cells, subcelldata);
+ GridTools::delete_unused_vertices (vertices, cells, subcelldata);
// ... and cells
GridReordering<dim>::invert_all_cells_of_negative_grid (vertices, cells);
GridReordering<dim>::reorder_cells (cells);
AssertThrow (in, ExcIO());
// do some clean-up on vertices...
- delete_unused_vertices (vertices, cells, subcelldata);
+ GridTools::delete_unused_vertices (vertices, cells, subcelldata);
// ...and cells
GridReordering<dim>::invert_all_cells_of_negative_grid (vertices, cells);
GridReordering<dim>::reorder_cells (cells);
AssertThrow (in, ExcIO());
// do some clean-up on vertices...
- delete_unused_vertices (vertices, cells, subcelldata);
+ GridTools::delete_unused_vertices (vertices, cells, subcelldata);
// ... and cells
GridReordering<2>::invert_all_cells_of_negative_grid (vertices, cells);
GridReordering<2>::reorder_cells (cells);
AssertThrow (in, ExcIO());
// do some clean-up on vertices...
- delete_unused_vertices (vertices, cells, subcelldata);
+ GridTools::delete_unused_vertices (vertices, cells, subcelldata);
// ... and cells
GridReordering<3>::invert_all_cells_of_negative_grid (vertices, cells);
GridReordering<3>::reorder_cells (cells);
AssertThrow (in, ExcIO());
// do some clean-up on vertices...
- delete_unused_vertices (vertices, cells, subcelldata);
+ GridTools::delete_unused_vertices (vertices, cells, subcelldata);
// ... and cells
GridReordering<dim>::invert_all_cells_of_negative_grid (vertices, cells);
GridReordering<dim>::reorder_cells (cells);
}
SubCellData subcelldata;
- delete_unused_vertices(vertices, cells, subcelldata);
+ GridTools::delete_unused_vertices(vertices, cells, subcelldata);
GridReordering<dim>::reorder_cells (cells);
tria->create_triangulation (vertices, cells, subcelldata);
#endif
cells[cell].vertices[i]=vertex_indices[cell*vertices_per_hex+i];
SubCellData subcelldata;
- delete_unused_vertices(vertices, cells, subcelldata);
+ GridTools::delete_unused_vertices(vertices, cells, subcelldata);
GridReordering<dim>::reorder_cells (cells);
tria->create_triangulation (vertices, cells, subcelldata);
#endif
-template <int dim>
-void
-GridIn<dim>::delete_unused_vertices (std::vector<Point<dim> > &vertices,
- std::vector<CellData<dim> > &cells,
- SubCellData &subcelldata)
-{
- // first check which vertices are
- // actually used
- std::vector<bool> vertex_used (vertices.size(), false);
- for (unsigned int c=0; c<cells.size(); ++c)
- for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
- vertex_used[cells[c].vertices[v]] = true;
-
- // then renumber the vertices that
- // are actually used in the same
- // order as they were beforehand
- const unsigned int invalid_vertex = deal_II_numbers::invalid_unsigned_int;
- std::vector<unsigned int> new_vertex_numbers (vertices.size(), invalid_vertex);
- unsigned int next_free_number = 0;
- for (unsigned int i=0; i<vertices.size(); ++i)
- if (vertex_used[i] == true)
- {
- new_vertex_numbers[i] = next_free_number;
- ++next_free_number;
- };
-
- // next replace old vertex numbers
- // by the new ones
- 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]];
-
- // same for boundary data
- for (unsigned int c=0; c<subcelldata.boundary_lines.size(); ++c)
- for (unsigned int v=0; v<GeometryInfo<1>::vertices_per_cell; ++v)
- subcelldata.boundary_lines[c].vertices[v]
- = new_vertex_numbers[subcelldata.boundary_lines[c].vertices[v]];
- for (unsigned int c=0; c<subcelldata.boundary_quads.size(); ++c)
- for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_cell; ++v)
- subcelldata.boundary_quads[c].vertices[v]
- = new_vertex_numbers[subcelldata.boundary_quads[c].vertices[v]];
-
- // finally copy over the vertices
- // which we really need to a new
- // array and replace the old one by
- // the new one
- std::vector<Point<dim> > tmp;
- tmp.reserve (std::count(vertex_used.begin(), vertex_used.end(), true));
- for (unsigned int v=0; v<vertices.size(); ++v)
- if (vertex_used[v] == true)
- tmp.push_back (vertices[v]);
- swap (vertices, tmp);
-}
-
-
-
template <int dim>
void GridIn<dim>::debug_output_grid (const std::vector<CellData<dim> > &/*cells*/,
const std::vector<Point<dim> > &/*vertices*/,
#endif
+template <int dim>
+void
+GridTools::delete_unused_vertices (std::vector<Point<dim> > &vertices,
+ std::vector<CellData<dim> > &cells,
+ SubCellData &subcelldata)
+{
+ // first check which vertices are
+ // actually used
+ std::vector<bool> vertex_used (vertices.size(), false);
+ for (unsigned int c=0; c<cells.size(); ++c)
+ for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+ vertex_used[cells[c].vertices[v]] = true;
+
+ // then renumber the vertices that
+ // are actually used in the same
+ // order as they were beforehand
+ const unsigned int invalid_vertex = deal_II_numbers::invalid_unsigned_int;
+ std::vector<unsigned int> new_vertex_numbers (vertices.size(), invalid_vertex);
+ unsigned int next_free_number = 0;
+ for (unsigned int i=0; i<vertices.size(); ++i)
+ if (vertex_used[i] == true)
+ {
+ new_vertex_numbers[i] = next_free_number;
+ ++next_free_number;
+ };
+
+ // next replace old vertex numbers
+ // by the new ones
+ 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]];
+
+ // same for boundary data
+ for (unsigned int c=0; c<subcelldata.boundary_lines.size(); ++c)
+ for (unsigned int v=0; v<GeometryInfo<1>::vertices_per_cell; ++v)
+ subcelldata.boundary_lines[c].vertices[v]
+ = new_vertex_numbers[subcelldata.boundary_lines[c].vertices[v]];
+ for (unsigned int c=0; c<subcelldata.boundary_quads.size(); ++c)
+ for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_cell; ++v)
+ subcelldata.boundary_quads[c].vertices[v]
+ = new_vertex_numbers[subcelldata.boundary_quads[c].vertices[v]];
+
+ // finally copy over the vertices
+ // which we really need to a new
+ // array and replace the old one by
+ // the new one
+ std::vector<Point<dim> > tmp;
+ tmp.reserve (std::count(vertex_used.begin(), vertex_used.end(), true));
+ for (unsigned int v=0; v<vertices.size(); ++v)
+ if (vertex_used[v] == true)
+ tmp.push_back (vertices[v]);
+ swap (vertices, tmp);
+}
+
+
// define some transformations in an anonymous namespace
#ifdef DEAL_II_ANON_NAMESPACE_BOGUS_WARNING
GridTools::diameter<deal_II_dimension> (const Triangulation<deal_II_dimension> &);
#endif
+template
+void GridTools::delete_unused_vertices (std::vector<Point<deal_II_dimension> > &,
+ std::vector<CellData<deal_II_dimension> > &,
+ SubCellData &);
+
template
void GridTools::shift<deal_II_dimension> (const Point<deal_II_dimension> &,
Triangulation<deal_II_dimension> &);
<h3>deal.II</h3>
<ol>
+ <li> <p>
+ New: There is now a new <code
+ class="member">GridTools::delete_unused_vertices</code>
+ function. Previously a <tt>private</tt> function in <code
+ class="class">GridIn</code> it has now been moved to and made
+ <tt>public</tt> in <code class="class">GridTools</code>.
+ <br>
+ (RH 2005/09/28)
+ </p>
+
<li> <p>
New: The <code
class="member">GridIn<dim>::read_netcdf(string