--- /dev/null
+/*---------------------------- grid_in.h ---------------------------*/
+/* $Id$ */
+#ifndef __grid_in_H
+#define __grid_in_H
+/*---------------------------- grid_in.h ---------------------------*/
+
+
+
+#include <base/exceptions.h>
+#include <base/smartpointer.h>
+#include <grid/forward_declarations.h>
+#include <iostream>
+#include <vector>
+#include <string>
+
+
+
+/**
+ * This class implements an input mechanism for grid data. It allows to
+ * read a grid structure into a triangulation object. Future versions
+ * will also allow to read data on this grid into vectors.
+ *
+ * At present, only UCD (unstructured cell data) is supported as input
+ * format for grid data. Any numerical data after the block of topological
+ * information is ignored.
+ *
+ * To read grid data, the triangulation to be fed with has to be empty.
+ * When giving a file which does not contain the assumed information or
+ * which does not keep to the right format, the state of the triangulation
+ * will be undefined afterwards. Upon input, only lines in one dimension
+ * and line and quads in two dimensions are accepted. All other cell types
+ * (e.g. triangles in two dimensions, quads and hexes in 3d) are rejected.
+ * The vertex and cell numbering in the UCD file, which
+ * need not be consecutively, is lost upon transfer to the triangulation
+ * object, since this one needs consecutively numbered elements.
+ *
+ * Material indicators are accepted to denote the material id of cells and
+ * to denote boundary part indication for lines in 2D. Read the according
+ * sections in the documentation of the \Ref{Triangulation} class for
+ * further details.
+ *
+ *
+ * \subsection{Structure of input grid data}
+ *
+ * It is your duty to use a correct numbering of vertices in the cell list,
+ * i.e. for lines, you have to first give the vertex with the lower coordinate
+ * value, then that with the higher coordinate value. For quadrilaterals in
+ * two dimensions, the vertex indices in the #quad# list have to be such that
+ * the vertices are numbered in counter-clockwise sense.
+ *
+ * In two dimensions, another difficulty occurs, which has to do with the sense
+ * of a quadrilateral. A quad consists of four lines which have a direction,
+ * which is per definitionem as follows:
+ * \begin{verbatim}
+ * 3-->--2
+ * | |
+ * ^ ^
+ * | |
+ * 0-->--1
+ * \end{verbatim}
+ * Now, two adjacent cells must have a vertex numbering such that the direction
+ * of the common side is the same. For example, the following two quads
+ * \begin{verbatim}
+ * 3---4---5
+ * | | |
+ * 0---1---2
+ * \end{verbatim}
+ * may be characterised by the vertex numbers (0 1 4 3) and (1 2 5 4), since
+ * the middle line would get the direction #1->4# when viewed from both cells.
+ * The numbering (0 1 4 3) and (5 4 1 2) would not be allowed, since the left
+ * quad would give the common line the direction #1->4#, while the right one
+ * would want to use #4->1#, leading to ambiguity. The #Triangulation# object
+ * is capable of detecting this special case, which can be eliminated by
+ * rotating the indices of the right quad by two. However, it would not
+ * know what to do if you gave the vertex indices (4 1 2 5), since then it
+ * would have to rotate by one element or three, the decision which to take is
+ * not yet implemented.
+ *
+ * There are more ambiguous cases, where the triangulation may not know what
+ * to do at all without the use of very sophisticated algorithms. On such example
+ * is the following:
+ * \begin{verbatim}
+ * 9---10-----11
+ * | | / |
+ * 6---7---8 |
+ * | | | |
+ * 3---4---5 |
+ * | | \ |
+ * 0---1------2
+ * \end{verbatim}
+ * Assume that you had numbered the vertices in the cells at the left boundary
+ * in a way, that the following line directions are induced:
+ * \begin{verbatim}
+ * 9->-10-----11
+ * ^ ^ / |
+ * 6->-7---8 |
+ * ^ ^ | |
+ * 3->-4---5 |
+ * ^ ^ \ |
+ * 0->-1------2
+ * \end{verbatim}
+ * (This could for example be done by using the indices (0 1 4 3), (3 4 7 6),
+ * (6 7 10 9) for the three cells). Now, you will not find a way of giving
+ * indices for the right cells, without introducing either ambiguity for
+ * one line or other, or without violating that within each cells, there must be
+ * one vertex from which both lines are directed away and the opposite one to
+ * which both adjacent lines point to.
+ *
+ * The solution in this case is to renumber one of the three left cells, e.g.
+ * by reverting the sense of the line between vertices 7 and 10 by numbering
+ * the top left cell by (9 6 7 10).
+ *
+ * But this is a thing that the triangulation
+ * object can't do for you, since it would involve backtracking to cells
+ * already created when we find that we can't number the indices of one of
+ * the rightmost cells consistently. It is neither clear how to do this
+ * backtracking nor whether it can be done with a stopping algorithm, if
+ * possible within polynomial time. This kind of numbering must be made
+ * upon construction of the coarse grid, unfortunately.
+ *
+ * @author Wolfgang Bangerth, 1998
+ */
+template <int dim>
+class GridIn
+{
+ public:
+ /**
+ * Constructor.
+ */
+ GridIn ();
+
+ /**
+ * Attach this triangulation
+ * to be fed with the grid data.
+ */
+ void attach_triangulation (Triangulation<dim> &tria);
+
+ /**
+ * Read grid data from an ucd file.
+ * Numerical data is ignored.
+ */
+ void read_ucd (istream &);
+
+ /**
+ * Exception
+ */
+ DeclException1 (ExcUnknownIdentifier,
+ string,
+ << "The identifier <" << arg1 << "> as name of a "
+ << "part in an UCD input file is unknown or the "
+ << "respective input routine is not implemented.");
+ /**
+ * Exception
+ */
+ DeclException0 (ExcNoTriangulationSelected);
+ /**
+ * Exception
+ */
+ DeclException2 (ExcInvalidVertexIndex,
+ int, int,
+ << "Trying to access invalid vertex index " << arg2
+ << " while creating cell " << arg1);
+ /**
+ * Exception
+ */
+ DeclException0 (ExcInternalError);
+ /**
+ * Exception
+ */
+ DeclException0 (ExcIO);
+
+ private:
+ /**
+ * Store address of the triangulation to
+ * be fed with the data read in.
+ */
+ SmartPointer<Triangulation<dim> > tria;
+};
+
+
+
+
+
+
+/*---------------------------- grid_in.h ---------------------------*/
+/* end of #ifndef __grid_in_H */
+#endif
+/*---------------------------- grid_in.h ---------------------------*/
-/**
- * This class implements an input mechanism for grid data. It allows to
- * read a grid structure into a triangulation object. Future versions
- * will also allow to read data on this grid into vectors.
- *
- * At present, only UCD (unstructured cell data) is supported as input
- * format for grid data. Any numerical data after the block of topological
- * information is ignored.
- *
- * To read grid data, the triangulation to be fed with has to be empty.
- * When giving a file which does not contain the assumed information or
- * which does not keep to the right format, the state of the triangulation
- * will be undefined afterwards. Upon input, only lines in one dimension
- * and line and quads in two dimensions are accepted. All other cell types
- * (e.g. triangles in two dimensions, quads and hexes in 3d) are rejected.
- * The vertex and cell numbering in the UCD file, which
- * need not be consecutively, is lost upon transfer to the triangulation
- * object, since this one needs consecutively numbered elements.
- *
- * Material indicators are accepted to denote the material id of cells and
- * to denote boundary part indication for lines in 2D. Read the according
- * sections in the documentation of the \Ref{Triangulation} class for
- * further details.
- *
- *
- * \subsection{Structure of input grid data}
- *
- * It is your duty to use a correct numbering of vertices in the cell list,
- * i.e. for lines, you have to first give the vertex with the lower coordinate
- * value, then that with the higher coordinate value. For quadrilaterals in
- * two dimensions, the vertex indices in the #quad# list have to be such that
- * the vertices are numbered in counter-clockwise sense.
- *
- * In two dimensions, another difficulty occurs, which has to do with the sense
- * of a quadrilateral. A quad consists of four lines which have a direction,
- * which is per definitionem as follows:
- * \begin{verbatim}
- * 3-->--2
- * | |
- * ^ ^
- * | |
- * 0-->--1
- * \end{verbatim}
- * Now, two adjacent cells must have a vertex numbering such that the direction
- * of the common side is the same. For example, the following two quads
- * \begin{verbatim}
- * 3---4---5
- * | | |
- * 0---1---2
- * \end{verbatim}
- * may be characterised by the vertex numbers (0 1 4 3) and (1 2 5 4), since
- * the middle line would get the direction #1->4# when viewed from both cells.
- * The numbering (0 1 4 3) and (5 4 1 2) would not be allowed, since the left
- * quad would give the common line the direction #1->4#, while the right one
- * would want to use #4->1#, leading to ambiguity. The #Triangulation# object
- * is capable of detecting this special case, which can be eliminated by
- * rotating the indices of the right quad by two. However, it would not
- * know what to do if you gave the vertex indices (4 1 2 5), since then it
- * would have to rotate by one element or three, the decision which to take is
- * not yet implemented.
- *
- * There are more ambiguous cases, where the triangulation may not know what
- * to do at all without the use of very sophisticated algorithms. On such example
- * is the following:
- * \begin{verbatim}
- * 9---10-----11
- * | | / |
- * 6---7---8 |
- * | | | |
- * 3---4---5 |
- * | | \ |
- * 0---1------2
- * \end{verbatim}
- * Assume that you had numbered the vertices in the cells at the left boundary
- * in a way, that the following line directions are induced:
- * \begin{verbatim}
- * 9->-10-----11
- * ^ ^ / |
- * 6->-7---8 |
- * ^ ^ | |
- * 3->-4---5 |
- * ^ ^ \ |
- * 0->-1------2
- * \end{verbatim}
- * (This could for example be done by using the indices (0 1 4 3), (3 4 7 6),
- * (6 7 10 9) for the three cells). Now, you will not find a way of giving
- * indices for the right cells, without introducing either ambiguity for
- * one line or other, or without violating that within each cells, there must be
- * one vertex from which both lines are directed away and the opposite one to
- * which both adjacent lines point to.
- *
- * The solution in this case is to renumber one of the three left cells, e.g.
- * by reverting the sense of the line between vertices 7 and 10 by numbering
- * the top left cell by (9 6 7 10).
- *
- * But this is a thing that the triangulation
- * object can't do for you, since it would involve backtracking to cells
- * already created when we find that we can't number the indices of one of
- * the rightmost cells consistently. It is neither clear how to do this
- * backtracking nor whether it can be done with a stopping algorithm, if
- * possible within polynomial time. This kind of numbering must be made
- * upon construction of the coarse grid, unfortunately.
- *
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim>
-class DataIn {
- public:
- /**
- * Constructor.
- */
- DataIn ();
-
- /**
- * Attach this triangulation
- * to be fed with the grid data.
- */
- void attach_triangulation (Triangulation<dim> *tria);
-
- /**
- * Read grid data from an ucd file.
- * Numerical data is ignored.
- */
- void read_ucd (istream &);
-
- /**
- * Exception
- */
- DeclException1 (ExcUnknownIdentifier,
- string,
- << "The identifier <" << arg1 << "> as name of a "
- << "part in an UCD input file is unknown or the "
- << "respective input routine is not implemented.");
- /**
- * Exception
- */
- DeclException0 (ExcNoTriangulationSelected);
- /**
- * Exception
- */
- DeclException2 (ExcInvalidVertexIndex,
- int, int,
- << "Trying to access invalid vertex index " << arg2
- << " while creating cell " << arg1);
- /**
- * Exception
- */
- DeclException0 (ExcInternalError);
- /**
- * Exception
- */
- DeclException0 (ExcIO);
-
- private:
- /**
- * Store address of the triangulation to
- * be fed with the data read in.
- */
- Triangulation<dim> *tria;
-};
-
-
-
-
-
-
-
/**
* This class is deprecated. Use the #DataOut# class instead.
*
--- /dev/null
+/* $Id$ */
+
+#include <grid/grid_in.h>
+#include <grid/tria.h>
+
+#include <map>
+
+
+
+template <int dim>
+GridIn<dim>::GridIn () :
+ tria(0) {};
+
+
+
+template <int dim>
+void GridIn<dim>::attach_triangulation (Triangulation<dim> &t)
+{
+ tria = &t;
+};
+
+
+
+template <int dim>
+void GridIn<dim>::read_ucd (istream &in)
+{
+ Assert (tria != 0, ExcNoTriangulationSelected());
+ Assert ((1<=dim) && (dim<=2), ExcNotImplemented());
+
+ AssertThrow (in, ExcIO());
+
+ // skip comments at start of file
+ char c;
+ while (in.get(c), c=='#')
+ {
+ char line[256];
+ in.get (line, 255, '\n'); // ignore rest of line, at most 256 chars
+ in.get (c); // ignore '\n' at end of line.
+ };
+
+ // put back first character of
+ // first non-comment line
+ in.putback (c);
+
+
+ unsigned int n_vertices;
+ unsigned int n_cells;
+ int dummy;
+
+ in >> n_vertices
+ >> n_cells
+ >> dummy // number of data vectors
+ >> dummy // cell data
+ >> dummy; // model data
+
+ // set up array of vertices
+ vector<Point<dim> > vertices (n_vertices);
+ // set up mapping between numbering
+ // in ucd-file (key) and in the
+ // vertices vector
+ map<int,int> vertex_indices;
+
+ for (unsigned int vertex=0; vertex<n_vertices; ++vertex)
+ {
+ int vertex_number;
+ double x[3];
+
+ // read vertex
+ in >> vertex_number
+ >> x[0] >> x[1] >> x[2];
+
+ // store vertex
+ for (unsigned int d=0; d<dim; ++d)
+ vertices[vertex](d) = x[d];
+ // store mapping; note that
+ // vertices_indices[i] is automatically
+ // created upon first usage
+ vertex_indices[vertex_number] = vertex;
+ };
+
+ // set up array of cells
+ vector<CellData<dim> > cells;
+ SubCellData subcelldata;
+
+ for (unsigned int cell=0; cell<n_cells; ++cell)
+ {
+ string cell_type;
+ int material_id;
+
+ in >> dummy // cell number
+ >> material_id;
+ in >> cell_type;
+
+ if (((cell_type == "line") && (dim == 1)) ||
+ ((cell_type == "quad") && (dim == 2)))
+ // found a cell
+ {
+ // allocate and read indices
+ cells.push_back (CellData<dim>());
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ in >> cells.back().vertices[i];
+ cells.back().material_id = material_id;
+
+ // transform from ucd to
+ // consecutive numbering
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ if (vertex_indices.find (cells.back().vertices[i]) != vertex_indices.end())
+ // vertex with this index exists
+ cells.back().vertices[i] = vertex_indices[cells.back().vertices[i]];
+ else
+ {
+ // no such vertex index
+ Assert (false, ExcInvalidVertexIndex(cell, cells.back().vertices[i]));
+ cells.back().vertices[i] = -1;
+ };
+ }
+ else
+ if ((cell_type == "line") && (dim == 2))
+ // boundary info
+ {
+ subcelldata.boundary_lines.push_back (CellData<1>());
+ in >> subcelldata.boundary_lines.back().vertices[0]
+ >> subcelldata.boundary_lines.back().vertices[1];
+ subcelldata.boundary_lines.back().material_id = material_id;
+
+ // transform from ucd to
+ // consecutive numbering
+ for (unsigned int i=0; i<2; ++i)
+ if (vertex_indices.find (subcelldata.boundary_lines.back().vertices[i]) !=
+ vertex_indices.end())
+ // vertex with this index exists
+ subcelldata.boundary_lines.back().vertices[i]
+ = vertex_indices[subcelldata.boundary_lines.back().vertices[i]];
+ else
+ {
+ // no such vertex index
+ Assert (false,
+ ExcInvalidVertexIndex(cell,
+ subcelldata.boundary_lines.back().vertices[i]));
+ subcelldata.boundary_lines.back().vertices[i] = -1;
+ };
+ }
+ else
+ // cannot read this
+ Assert (false, ExcUnknownIdentifier(cell_type));
+ };
+
+ // check that no forbidden arrays are used
+ Assert (subcelldata.check_consistency(dim), ExcInternalError());
+
+ AssertThrow (in, ExcIO());
+
+ tria->create_triangulation (vertices, cells, subcelldata);
+};
+
+
+
+//explicit instantiations
+template class GridIn<deal_II_dimension>;
-template <int dim>
-DataIn<dim>::DataIn () :
- tria(0) {};
-
-
-
-template <int dim>
-void DataIn<dim>::attach_triangulation (Triangulation<dim> *t) {
- tria = t;
-};
-
-
-
-template <int dim>
-void DataIn<dim>::read_ucd (istream &in) {
- Assert (tria != 0, ExcNoTriangulationSelected());
- Assert ((1<=dim) && (dim<=2), ExcNotImplemented());
-
- AssertThrow (in, ExcIO());
-
- // skip comments at start of file
- char c;
- while (in.get(c), c=='#')
- {
- char line[256];
- in.get (line, 255, '\n'); // ignore rest of line, at most 256 chars
- in.get (c); // ignore '\n' at end of line.
- };
-
- // put back first character of
- // first non-comment line
- in.putback (c);
-
-
- unsigned int n_vertices;
- unsigned int n_cells;
- int dummy;
-
- in >> n_vertices
- >> n_cells
- >> dummy // number of data vectors
- >> dummy // cell data
- >> dummy; // model data
-
- // set up array of vertices
- vector<Point<dim> > vertices (n_vertices);
- // set up mapping between numbering
- // in ucd-file (key) and in the
- // vertices vector
- map<int,int> vertex_indices;
-
- for (unsigned int vertex=0; vertex<n_vertices; ++vertex)
- {
- int vertex_number;
- double x[3];
-
- // read vertex
- in >> vertex_number
- >> x[0] >> x[1] >> x[2];
-
- // store vertex
- for (unsigned int d=0; d<dim; ++d)
- vertices[vertex](d) = x[d];
- // store mapping; note that
- // vertices_indices[i] is automatically
- // created upon first usage
- vertex_indices[vertex_number] = vertex;
- };
-
- // set up array of cells
- vector<CellData<dim> > cells;
- SubCellData subcelldata;
-
- for (unsigned int cell=0; cell<n_cells; ++cell)
- {
- string cell_type;
- int material_id;
-
- in >> dummy // cell number
- >> material_id;
- in >> cell_type;
-
- if (((cell_type == "line") && (dim == 1)) ||
- ((cell_type == "quad") && (dim == 2)))
- // found a cell
- {
- // allocate and read indices
- cells.push_back (CellData<dim>());
- for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
- in >> cells.back().vertices[i];
- cells.back().material_id = material_id;
-
- // transform from ucd to
- // consecutive numbering
- for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
- if (vertex_indices.find (cells.back().vertices[i]) != vertex_indices.end())
- // vertex with this index exists
- cells.back().vertices[i] = vertex_indices[cells.back().vertices[i]];
- else
- {
- // no such vertex index
- Assert (false, ExcInvalidVertexIndex(cell, cells.back().vertices[i]));
- cells.back().vertices[i] = -1;
- };
- }
- else
- if ((cell_type == "line") && (dim == 2))
- // boundary info
- {
- subcelldata.boundary_lines.push_back (CellData<1>());
- in >> subcelldata.boundary_lines.back().vertices[0]
- >> subcelldata.boundary_lines.back().vertices[1];
- subcelldata.boundary_lines.back().material_id = material_id;
-
- // transform from ucd to
- // consecutive numbering
- for (unsigned int i=0; i<2; ++i)
- if (vertex_indices.find (subcelldata.boundary_lines.back().vertices[i]) !=
- vertex_indices.end())
- // vertex with this index exists
- subcelldata.boundary_lines.back().vertices[i]
- = vertex_indices[subcelldata.boundary_lines.back().vertices[i]];
- else
- {
- // no such vertex index
- Assert (false,
- ExcInvalidVertexIndex(cell,
- subcelldata.boundary_lines.back().vertices[i]));
- subcelldata.boundary_lines.back().vertices[i] = -1;
- };
- }
- else
- // cannot read this
- Assert (false, ExcUnknownIdentifier(cell_type));
- };
-
- // check that no forbidden arrays are used
- Assert (subcelldata.check_consistency(dim), ExcInternalError());
-
- AssertThrow (in, ExcIO());
-
- tria->create_triangulation (vertices, cells, subcelldata);
-};
-
-
-
-
-
-
template <int dim>
DataOut_Old<dim>::DataEntry::DataEntry () :
//explicit instantiations
-
-template class DataIn<deal_II_dimension>;
template class DataOut_Old<deal_II_dimension>;
kdoc.library-files.rerun = $(kdoc.library-files:%.kdoc=%.kdoc.rerun)
-KDOCFLAGS = -I../scripts/kdoc ../scripts/kdoc/kdoc -a -p
-#KDOCFLAGS = /home/people/wolf/tmp/kdoc/kdoc-2.a22/bin/kdoc -a -p
+#KDOCFLAGS = -I../scripts/kdoc ../scripts/kdoc/kdoc -a -p
+KDOCFLAGS = /home/people/wolf/tmp/kdoc/kdoc-2.a22/bin/kdoc -a -p
ifeq ($(shell uname),SunOS)
DOCPP = /home/people/wolf/Config/doc++/doc++
#DOCPP = /home/people/kanschat/bin/doc++
This is possible: <acronym>deal.II</acronym> offers the possibility of
reading a complete triangulation from a file in the <tt>ucd</tt>-format
used by avs. A <tt>ucd</tt>-file can be read with the function<br>
-<code>void DataIn::read_ucd(istream&)</code><br>
+<code>void GridIn::read_ucd(istream&)</code><br>
</p>
<p><span class="parhead">
may have the wrong sign (plus instead of minus or vice versa).
A more detailed description of the problems
encountered in two dimensions can be found in the
-<a href="http://gaia.iwr.uni-heidelberg.de/~deal/doc/auto/kdoc/basic/DataIn.html" target="_top"> <code>DataIn</code> class description</a>.
+<a href="../../auto/kdoc/grid/GridIn.html" target="_top"> <code>GridIn</code> class description</a>.
</p>