/**
* Read grid data from a unstructured vtk file. The vtk file may contain
- * the following VTK cell types: VTK_HEXAHEDRON, VTK_QUAD, and VTK_LINE.
+ * the following VTK cell types: VTK_HEXAHEDRON (12), VTK_TETRA (10),
+ * VTK_QUAD (9), VTK_TRIANGLE (5), and VTK_LINE (3).
*
* Depending on the template dimension, only some of the above are accepted.
*
* In particular, in three dimensions, this function expects the file to
* contain
*
- * - VTK_HEXAHEDRON cell types
- * - VTK_QUAD cell types, to specify optional boundary or interior quad faces
+ * - VTK_HEXAHEDRON/VTK_TETRA cell types
+ * - VTK_QUAD/VTK_TRIANGLE cell types, to specify optional boundary or
+ * interior quad faces
* - VTK_LINE cell types, to specify optional boundary or interior edges
*
* In two dimensions:
*
- * - VTK_QUAD cell types
+ * - VTK_QUAD/VTK_TRIANGLE cell types
* - VTK_LINE cell types, to specify optional boundary or interior edges
*
* In one dimension
*
* The companion GridOut::write_vtk function can be used to write VTK files
* compatible with this method.
+ *
+ * @ingroup simplex
*/
void
read_vtk(std::istream &in);
unsigned int n_geometric_objects = 0;
unsigned int n_ints;
+ bool is_quad_or_hex_mesh = false;
+ bool is_tria_or_tet_mesh = false;
+
if (keyword == "CELLS")
{
+ // jump to the `CELL_TYPES` section and read in cell types
+ std::vector<unsigned int> cell_types;
+ {
+ std::streampos oldpos = in.tellg();
+
+
+ while (in >> keyword)
+ if (keyword == "CELL_TYPES")
+ {
+ in >> n_ints;
+
+ cell_types.resize(n_ints);
+
+ for (unsigned int i = 0; i < n_ints; ++i)
+ in >> cell_types[i];
+
+ break;
+ }
+
+ in.seekg(oldpos);
+ }
+
in >> n_geometric_objects;
in >> n_ints; // Ignore this, since we don't need it.
unsigned int type;
in >> type;
- if (type == 8)
+ if (cell_types[count] == 10 || cell_types[count] == 12)
{
+ if (cell_types[count] == 10)
+ is_tria_or_tet_mesh = true;
+ if (cell_types[count] == 12)
+ is_quad_or_hex_mesh = true;
+
// we assume that the file contains first all cells,
// and only then any faces or lines
AssertThrow(subcelldata.boundary_quads.size() == 0 &&
subcelldata.boundary_lines.size() == 0,
ExcNotImplemented());
- cells.emplace_back();
+ cells.emplace_back(type);
for (unsigned int j = 0; j < type; j++) // loop to feed data
in >> cells.back().vertices[j];
cells.back().material_id = 0;
}
- else if (type == 4)
+ else if (cell_types[count] == 5 || cell_types[count] == 9)
{
+ if (cell_types[count] == 5)
+ is_tria_or_tet_mesh = true;
+ if (cell_types[count] == 9)
+ is_quad_or_hex_mesh = true;
+
// we assume that the file contains first all cells,
// then all faces, and finally all lines
AssertThrow(subcelldata.boundary_lines.size() == 0,
ExcNotImplemented());
- subcelldata.boundary_quads.emplace_back();
+ subcelldata.boundary_quads.emplace_back(type);
for (unsigned int j = 0; j < type;
j++) // loop to feed the data to the boundary
subcelldata.boundary_quads.back().material_id = 0;
}
- else if (type == 2)
+ else if (cell_types[count] == 3)
{
- subcelldata.boundary_lines.emplace_back();
+ subcelldata.boundary_lines.emplace_back(type);
for (unsigned int j = 0; j < type;
j++) // loop to feed the data to the boundary
"While reading VTK file, unknown file type encountered"));
}
}
-
else if (dim == 2)
{
for (unsigned int count = 0; count < n_geometric_objects; count++)
unsigned int type;
in >> type;
- if (type == 4)
+ if (cell_types[count] == 5 || cell_types[count] == 9)
{
// we assume that the file contains first all cells,
// and only then any faces
AssertThrow(subcelldata.boundary_lines.size() == 0,
ExcNotImplemented());
- cells.emplace_back();
+ if (cell_types[count] == 5)
+ is_tria_or_tet_mesh = true;
+ if (cell_types[count] == 9)
+ is_quad_or_hex_mesh = true;
+
+ cells.emplace_back(type);
for (unsigned int j = 0; j < type; j++) // loop to feed data
in >> cells.back().vertices[j];
cells.back().material_id = 0;
}
- else if (type == 2)
+ else if (cell_types[count] == 3)
{
// If this is encountered, the pointer comes out of the loop
// and starts processing boundaries.
- subcelldata.boundary_lines.emplace_back();
+ subcelldata.boundary_lines.emplace_back(type);
for (unsigned int j = 0; j < type;
j++) // loop to feed the data to the boundary
in >> type;
AssertThrow(
- type == 2,
+ cell_types[count] == 3 && type == 2,
ExcMessage(
"While reading VTK file, unknown cell type encountered"));
- cells.emplace_back();
+ cells.emplace_back(type);
for (unsigned int j = 0; j < type; j++) // loop to feed data
in >> cells.back().vertices[j];
Assert(subcelldata.check_consistency(dim), ExcInternalError());
- GridTools::delete_unused_vertices(vertices, cells, subcelldata);
- if (dim == spacedim)
- GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(
- vertices, cells);
+ // make sure that only either simplex or hypercube cells are available
+ //
+ // TODO: the functions below (GridTools::delete_unused_vertices(),
+ // GridTools::invert_all_cells_of_negative_grid(),
+ // GridReordering::reorder_cells(),
+ // Triangulation::create_triangulation_compatibility()) need to be
+ // revisited for simplex meshes
+ AssertThrow(dim == 1 || (is_tria_or_tet_mesh ^ is_quad_or_hex_mesh),
+ ExcNotImplemented());
+
+ if (dim == 1 || is_quad_or_hex_mesh)
+ {
+ GridTools::delete_unused_vertices(vertices, cells, subcelldata);
+
+ if (dim == spacedim)
+ GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(
+ vertices, cells);
- GridReordering<dim, spacedim>::reorder_cells(cells);
- tria->create_triangulation_compatibility(vertices, cells, subcelldata);
+ GridReordering<dim, spacedim>::reorder_cells(cells);
+ tria->create_triangulation_compatibility(vertices,
+ cells,
+ subcelldata);
- return;
+ return;
+ }
+ else
+ {
+ tria->create_triangulation(vertices, cells, subcelldata);
+ }
}
else
AssertThrow(false,
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Read a file in the MSH format with (linear) triangular and tetrahedral
+// elements created by the GMSH program.
+
+#include <deal.II/grid/grid_in.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria.h>
+
+#include <fstream>
+
+#include "../tests.h"
+
+template <int dim, int spacedim>
+void
+check_file(const std::string &file_name)
+{
+ Triangulation<dim, spacedim> tria;
+
+ GridIn<dim, spacedim> grid_in;
+ grid_in.attach_triangulation(tria);
+ std::ifstream input_file(file_name);
+ grid_in.read_vtk(input_file);
+
+ GridOut grid_out;
+#if false
+ std::ofstream out("mesh.out.vtk");
+ grid_out.write_vtk(tria, out);
+#else
+ grid_out.write_vtk(tria, deallog.get_file_stream());
+#endif
+
+ deallog << "OK!" << std::endl;
+}
+
+
+int
+main()
+{
+ initlog();
+ // TRIANGULAR ELEMENTS
+ // dim = spacedim = 2
+ deallog.push("triangluar_elements_dim2_spacedim2: ");
+ check_file<2, 2>(std::string(SOURCE_DIR "/grid_in_vtk/tri.vtk"));
+ deallog.pop();
+
+ // dim = 2, spacedim = 3
+ deallog.push("triangluar_elements_dim2_spacedim3: ");
+ check_file<2, 3>(std::string(SOURCE_DIR "/grid_in_vtk/tri.vtk"));
+ deallog.pop();
+
+
+
+ // TETRAHEDRAL ELEMENTS
+ // dim = spacedim = 3
+ deallog.push("tetrahedral_elements_dim3_spacedim3: ");
+ check_file<3, 3>(std::string(SOURCE_DIR "/grid_in_vtk/tet.vtk"));
+ deallog.pop();
+}
--- /dev/null
+
+# vtk DataFile Version 3.0
+Triangulation generated with deal.II
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 5 double
+0.00000 0.00000 0
+1.00000 0.00000 0
+1.00000 1.00000 0
+0.00000 1.00000 0
+0.500000 0.500000 0
+
+CELLS 4 16
+3 0 1 4
+3 3 0 4
+3 1 2 4
+3 2 3 4
+
+CELL_TYPES 4
+5 5 5 5
+
+
+
+CELL_DATA 4
+SCALARS MaterialID int 1
+LOOKUP_TABLE default
+0 0 0 0
+
+
+
+SCALARS ManifoldID int 1
+LOOKUP_TABLE default
+-1 -1 -1 -1
+
+
+DEAL:triangluar_elements_dim2_spacedim2: ::OK!
+# vtk DataFile Version 3.0
+Triangulation generated with deal.II
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 5 double
+0.00000 0.00000 0.00000
+1.00000 0.00000 0.00000
+1.00000 1.00000 0.00000
+0.00000 1.00000 0.00000
+0.500000 0.500000 0.00000
+
+CELLS 4 16
+3 0 1 4
+3 3 0 4
+3 1 2 4
+3 2 3 4
+
+CELL_TYPES 4
+5 5 5 5
+
+
+
+CELL_DATA 4
+SCALARS MaterialID int 1
+LOOKUP_TABLE default
+0 0 0 0
+
+
+
+SCALARS ManifoldID int 1
+LOOKUP_TABLE default
+-1 -1 -1 -1
+
+
+DEAL:triangluar_elements_dim2_spacedim3: ::OK!
+# vtk DataFile Version 3.0
+Triangulation generated with deal.II
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 14 double
+0.00000 0.00000 1.00000
+0.00000 0.00000 0.00000
+0.00000 1.00000 1.00000
+0.00000 1.00000 0.00000
+1.00000 0.00000 1.00000
+1.00000 0.00000 0.00000
+1.00000 1.00000 1.00000
+1.00000 1.00000 0.00000
+0.00000 0.500000 0.500000
+1.00000 0.500000 0.500000
+0.500000 0.00000 0.500000
+0.500000 1.00000 0.500000
+0.500000 0.500000 0.00000
+0.500000 0.500000 1.00000
+
+CELLS 24 120
+4 13 10 11 9
+4 9 10 11 12
+4 11 8 10 13
+4 8 11 10 12
+4 4 13 6 9
+4 3 8 11 2
+4 1 8 0 10
+4 0 8 2 13
+4 4 10 0 13
+4 9 4 10 5
+4 11 2 13 6
+4 9 11 6 7
+4 1 3 8 12
+4 11 3 7 12
+4 12 7 9 5
+4 12 10 1 5
+4 13 4 10 9
+4 8 11 2 13
+4 6 13 11 9
+4 8 0 10 13
+4 10 1 8 12
+4 11 8 3 12
+4 7 9 11 12
+4 12 9 10 5
+
+CELL_TYPES 24
+10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10
+
+
+
+CELL_DATA 24
+SCALARS MaterialID int 1
+LOOKUP_TABLE default
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+
+
+
+SCALARS ManifoldID int 1
+LOOKUP_TABLE default
+-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1
+
+
+DEAL:tetrahedral_elements_dim3_spacedim3: ::OK!
projects / dealii.git / commitdiff