25#include <boost/algorithm/string.hpp>
26#include <boost/archive/binary_iarchive.hpp>
27#include <boost/io/ios_state.hpp>
28#include <boost/property_tree/ptree.hpp>
29#include <boost/property_tree/xml_parser.hpp>
30#include <boost/serialization/serialization.hpp>
32#ifdef DEAL_II_GMSH_WITH_API
43#ifdef DEAL_II_WITH_ASSIMP
44# include <assimp/Importer.hpp>
45# include <assimp/postprocess.h>
46# include <assimp/scene.h>
49#ifdef DEAL_II_TRILINOS_WITH_SEACAS
76 template <
int spacedim>
78 assign_1d_boundary_ids(
85 if (cell->face(face_no)->at_boundary())
87 if (cell->face(face_no)->vertex(0) == pair.
first)
89 cell->face(face_no)->set_boundary_id(pair.second);
95 template <
int dim,
int spacedim>
97 assign_1d_boundary_ids(
109 template <
int dim,
int spacedim>
117 const auto n_hypercube_vertices =
118 ReferenceCells::get_hypercube<dim>().n_vertices();
119 bool is_only_hypercube =
true;
121 if (cell.
vertices.size() != n_hypercube_vertices)
123 is_only_hypercube =
false;
131 if (is_only_hypercube)
136template <
int dim,
int spacedim>
138 : tria(nullptr, typeid(*this).name())
139 , default_format(ucd)
144template <
int dim,
int spacedim>
146 : tria(&t, typeid(*this).name())
147 , default_format(ucd)
152template <
int dim,
int spacedim>
161template <
int dim,
int spacedim>
173 text[0] =
"# vtk DataFile Version 3.0";
176 text[3] =
"DATASET UNSTRUCTURED_GRID";
178 for (
unsigned int i = 0; i < 4; ++i)
183 line.compare(text[i]) == 0,
186 "While reading VTK file, failed to find a header line with text <") +
193 std::vector<Point<spacedim>>
vertices;
194 std::vector<CellData<dim>> cells;
203 if (keyword ==
"POINTS")
205 unsigned int n_vertices;
210 for (
unsigned int vertex = 0; vertex < n_vertices; ++vertex)
214 in >> x[0] >> x[1] >> x[2];
217 for (
unsigned int d = 0; d < spacedim; ++d)
225 "While reading VTK file, failed to find POINTS section"));
229 unsigned int n_geometric_objects = 0;
232 if (keyword ==
"CELLS")
235 std::vector<unsigned int> cell_types;
237 std::streampos oldpos = in.tellg();
240 while (in >> keyword)
241 if (keyword ==
"CELL_TYPES")
245 cell_types.resize(n_ints);
247 for (
unsigned int i = 0; i < n_ints; ++i)
256 in >> n_geometric_objects;
261 for (
unsigned int count = 0; count < n_geometric_objects; ++count)
263 unsigned int n_vertices;
267 if (cell_types[count] == 10 || cell_types[count] == 12)
275 cells.emplace_back(n_vertices);
277 for (
unsigned int j = 0; j < n_vertices;
279 in >> cells.back().vertices[j];
283 if (cell_types[count] == 12)
285 std::swap(cells.back().vertices[2],
286 cells.back().vertices[3]);
287 std::swap(cells.back().vertices[6],
288 cells.back().vertices[7]);
291 cells.back().material_id = 0;
294 else if (cell_types[count] == 5 || cell_types[count] == 9)
303 for (
unsigned int j = 0; j < n_vertices;
310 else if (cell_types[count] == 3)
314 for (
unsigned int j = 0; j < n_vertices;
325 "While reading VTK file, unknown cell type encountered"));
330 for (
unsigned int count = 0; count < n_geometric_objects; ++count)
332 unsigned int n_vertices;
336 if (cell_types[count] == 5 || cell_types[count] == 9)
343 cells.emplace_back(n_vertices);
345 for (
unsigned int j = 0; j < n_vertices;
347 in >> cells.back().vertices[j];
351 if (cell_types[count] == 9)
355 std::swap(cells.back().vertices[2],
356 cells.back().vertices[3]);
359 cells.back().material_id = 0;
362 else if (cell_types[count] == 3)
368 for (
unsigned int j = 0; j < n_vertices;
381 "While reading VTK file, unknown cell type encountered"));
386 for (
unsigned int count = 0; count < n_geometric_objects; ++count)
392 cell_types[count] == 3 && type == 2,
394 "While reading VTK file, unknown cell type encountered"));
395 cells.emplace_back(type);
397 for (
unsigned int j = 0; j < type; ++j)
398 in >> cells.back().vertices[j];
400 cells.back().material_id = 0;
406 "While reading VTK file, failed to find CELLS section"));
413 keyword ==
"CELL_TYPES",
415 "While reading VTK file, missing CELL_TYPES section. Found <" +
416 keyword +
"> instead.")));
420 n_ints == n_geometric_objects,
421 ExcMessage(
"The VTK reader found a CELL_DATA statement "
422 "that lists a total of " +
424 " cell data objects, but this needs to "
425 "equal the number of cells (which is " +
427 ") plus the number of quads (" +
429 " in 3d or the number of lines (" +
434 for (
unsigned int i = 0; i < n_ints; ++i)
438 while (in >> keyword)
439 if (keyword ==
"CELL_DATA")
445 ExcMessage(
"The VTK reader found a CELL_DATA statement "
446 "that lists a total of " +
448 " cell data objects, but this needs to "
449 "equal the number of cells (which is " +
451 ") plus the number of quads (" +
454 " in 3d or the number of lines (" +
459 const std::vector<std::string> data_sets{
"MaterialID",
462 for (
unsigned int i = 0; i < data_sets.size(); ++i)
465 while (in >> keyword)
466 if (keyword ==
"SCALARS")
471 std::string field_name;
473 if (std::find(data_sets.begin(),
475 field_name) == data_sets.end())
487 std::getline(in, line);
490 std::min(
static_cast<std::size_t
>(3),
491 line.size() - 1)) ==
"int",
493 "While reading VTK file, material- and manifold IDs can only have type 'int'."));
497 keyword ==
"LOOKUP_TABLE",
499 "While reading VTK file, missing keyword 'LOOKUP_TABLE'."));
503 keyword ==
"default",
505 "While reading VTK file, missing keyword 'default'."));
512 for (
unsigned int i = 0; i < cells.size(); ++i)
516 if (field_name ==
"MaterialID")
517 cells[i].material_id =
519 else if (field_name ==
"ManifoldID")
520 cells[i].manifold_id =
532 if (field_name ==
"MaterialID")
533 boundary_quad.material_id =
535 else if (field_name ==
"ManifoldID")
536 boundary_quad.manifold_id =
545 if (field_name ==
"MaterialID")
546 boundary_line.material_id =
548 else if (field_name ==
"ManifoldID")
549 boundary_line.manifold_id =
561 if (field_name ==
"MaterialID")
562 boundary_line.material_id =
564 else if (field_name ==
"ManifoldID")
565 boundary_line.manifold_id =
575 apply_grid_fixup_functions(
vertices, cells, subcelldata);
581 "While reading VTK file, failed to find CELLS section"));
586template <
int dim,
int spacedim>
590 namespace pt = boost::property_tree;
592 pt::read_xml(in, tree);
593 auto section = tree.get_optional<std::string>(
"VTKFile.dealiiData");
597 "While reading a VTU file, failed to find dealiiData section. "
598 "Notice that we can only read grid files in .vtu format that "
599 "were created by the deal.II library, using a call to "
600 "GridOut::write_vtu(), where the flag "
601 "GridOutFlags::Vtu::serialize_triangulation is set to true."));
605 const auto string_archive =
607 std::istringstream in_stream(string_archive);
608 boost::archive::binary_iarchive ia(in_stream);
613template <
int dim,
int spacedim>
617 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
621 skip_comment_lines(in,
'#');
632 "Expected '-1' before and after a section."));
644 std::getline(in, line);
646 boost::algorithm::trim(line);
647 if (line.compare(
"-1") == 0)
658 std::vector<Point<spacedim>>
vertices;
677 in >> dummy >> dummy >> dummy;
680 in >> x[0] >> x[1] >> x[2];
684 for (
unsigned int d = 0; d < spacedim; ++d)
699 AssertThrow(tmp == 2412, ExcUnknownSectionType(tmp));
701 std::vector<CellData<dim>> cells;
728 in >> type >> dummy >> dummy >> dummy >> dummy;
730 AssertThrow((type == 11) || (type == 44) || (type == 94) || (type == 115),
731 ExcUnknownElementType(type));
733 if ((((type == 44) || (type == 94)) && (dim == 2)) ||
734 ((type == 115) && (dim == 3)))
736 const auto reference_cell = ReferenceCells::get_hypercube<dim>();
737 cells.emplace_back();
742 .vertices[reference_cell.unv_vertex_to_deal_vertex(v)];
744 cells.back().material_id = 0;
747 cells.back().vertices[v] =
vertex_indices[cells.back().vertices[v]];
749 cell_indices[object_index] = n_cells;
753 else if (((type == 11) && (dim == 2)) ||
754 ((type == 11) && (dim == 3)))
757 in >> dummy >> dummy >> dummy;
762 for (
unsigned int &vertex :
768 for (
unsigned int &vertex :
772 line_indices[object_index] = n_lines;
776 else if (((type == 44) || (type == 94)) && (dim == 3))
787 .vertices[reference_cell.unv_vertex_to_deal_vertex(v)];
791 for (
unsigned int &vertex :
795 quad_indices[object_index] = n_quads;
803 "> when running in dim=" +
822 AssertThrow((tmp == 2467) || (tmp == 2477), ExcUnknownSectionType(tmp));
838 in >> dummy >> dummy >> dummy >> dummy >> dummy >> dummy >>
849 std::getline(in, line);
852 "The line before the line containing an ID has too "
853 "many entries. This is not a valid UNV file."));
855 std::getline(in, line);
857 std::istringstream id_stream(line);
860 !id_stream.fail() && id_stream.eof(),
862 "The given UNV file contains a boundary or material id set to '" +
864 "', which cannot be parsed as a fixed-width integer, whereas "
865 "deal.II only supports integer boundary and material ids. To fix "
866 "this, ensure that all such ids are given integer values."));
869 id <=
long(std::numeric_limits<types::material_id>::max()),
870 ExcMessage(
"The provided integer id '" + std::to_string(
id) +
871 "' is not convertible to either types::material_id nor "
872 "types::boundary_id."));
874 const unsigned int n_lines =
875 (n_entities % 2 == 0) ? (n_entities / 2) : ((n_entities + 1) / 2);
877 for (
unsigned int line = 0; line < n_lines; ++line)
879 unsigned int n_fragments;
881 if (line == n_lines - 1)
882 n_fragments = (n_entities % 2 == 0) ? (2) : (1);
886 for (
unsigned int no_fragment = 0; no_fragment < n_fragments;
890 in >> dummy >> no >> dummy >> dummy;
892 if (cell_indices.count(no) > 0)
893 cells[cell_indices[no]].material_id = id;
895 if (line_indices.count(no) > 0)
899 if (quad_indices.count(no) > 0)
907 apply_grid_fixup_functions(
vertices, cells, subcelldata);
913template <
int dim,
int spacedim>
916 const bool apply_all_indicators_to_manifolds)
918 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
922 skip_comment_lines(in,
'#');
925 unsigned int n_vertices;
926 unsigned int n_cells;
929 in >> n_vertices >> n_cells >> dummy
935 std::vector<Point<spacedim>>
vertices(n_vertices);
941 for (
unsigned int vertex = 0; vertex < n_vertices; ++vertex)
948 in >> vertex_number >> x[0] >> x[1] >> x[2];
951 for (
unsigned int d = 0; d < spacedim; ++d)
960 std::vector<CellData<dim>> cells;
963 for (
unsigned int cell = 0; cell < n_cells; ++cell)
972 std::string cell_type;
976 unsigned int material_id;
982 if (((cell_type ==
"line") && (dim == 1)) ||
983 ((cell_type ==
"quad") && (dim == 2)) ||
984 ((cell_type ==
"hex") && (dim == 3)))
988 cells.emplace_back();
993 Assert(material_id <= std::numeric_limits<types::material_id>::max(),
996 std::numeric_limits<types::material_id>::max()));
1001 if (apply_all_indicators_to_manifolds)
1002 cells.back().manifold_id =
1004 cells.back().material_id = material_id;
1012 cells.back().vertices[i] =
1018 ExcInvalidVertexIndex(cell,
1019 cells.back().vertices[i]));
1024 else if ((cell_type ==
"line") && ((dim == 2) || (dim == 3)))
1032 Assert(material_id <= std::numeric_limits<types::boundary_id>::max(),
1035 std::numeric_limits<types::boundary_id>::max()));
1044 if (apply_all_indicators_to_manifolds)
1061 for (
unsigned int &vertex :
1069 AssertThrow(
false, ExcInvalidVertexIndex(cell, vertex));
1073 else if ((cell_type ==
"quad") && (dim == 3))
1082 Assert(material_id <= std::numeric_limits<types::boundary_id>::max(),
1085 std::numeric_limits<types::boundary_id>::max()));
1094 if (apply_all_indicators_to_manifolds)
1111 for (
unsigned int &vertex :
1119 Assert(
false, ExcInvalidVertexIndex(cell, vertex));
1125 AssertThrow(
false, ExcUnknownIdentifier(cell_type));
1130 apply_grid_fixup_functions(
vertices, cells, subcelldata);
1136 template <
int dim,
int spacedim>
1143 read_in_abaqus(std::istream &in);
1145 write_out_avs_ucd(std::ostream &out)
const;
1148 const double tolerance;
1151 get_global_node_numbers(
const int face_cell_no,
1152 const int face_cell_face_no)
const;
1155 std::vector<std::vector<double>> node_list;
1158 std::vector<std::vector<double>> cell_list;
1160 std::vector<std::vector<double>> face_list;
1163 std::map<std::string, std::vector<int>> elsets_list;
1167template <
int dim,
int spacedim>
1170 const bool apply_all_indicators_to_manifolds)
1172 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
1177 Assert((spacedim == 2 && dim == spacedim) ||
1178 (spacedim == 3 && (dim == spacedim || dim == spacedim - 1)),
1184 Abaqus_to_UCD<dim, spacedim> abaqus_to_ucd;
1185 abaqus_to_ucd.read_in_abaqus(in);
1187 std::stringstream in_ucd;
1188 abaqus_to_ucd.write_out_avs_ucd(in_ucd);
1196 read_ucd(in_ucd, apply_all_indicators_to_manifolds);
1198 catch (std::exception &exc)
1200 std::cerr <<
"Exception on processing internal UCD data: " << std::endl
1201 << exc.what() << std::endl;
1206 "Internal conversion from ABAQUS file to UCD format was unsuccessful. "
1207 "More information is provided in an error message printed above. "
1208 "Are you sure that your ABAQUS mesh file conforms with the requirements "
1209 "listed in the documentation?"));
1216 "Internal conversion from ABAQUS file to UCD format was unsuccessful. "
1217 "Are you sure that your ABAQUS mesh file conforms with the requirements "
1218 "listed in the documentation?"));
1223template <
int dim,
int spacedim>
1227 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
1233 skip_comment_lines(in,
'#');
1239 AssertThrow(line ==
"MeshVersionFormatted 0", ExcInvalidDBMESHInput(line));
1241 skip_empty_lines(in);
1245 AssertThrow(line ==
"Dimension", ExcInvalidDBMESHInput(line));
1246 unsigned int dimension;
1248 AssertThrow(dimension == dim, ExcDBMESHWrongDimension(dimension));
1249 skip_empty_lines(in);
1262 while (getline(in, line), line.find(
"# END") == std::string::npos)
1264 skip_empty_lines(in);
1269 AssertThrow(line ==
"Vertices", ExcInvalidDBMESHInput(line));
1271 unsigned int n_vertices;
1275 std::vector<Point<spacedim>>
vertices(n_vertices);
1276 for (
unsigned int vertex = 0; vertex < n_vertices; ++vertex)
1279 for (
unsigned int d = 0; d < dim; ++d)
1286 skip_empty_lines(in);
1292 AssertThrow(line ==
"Edges", ExcInvalidDBMESHInput(line));
1294 unsigned int n_edges;
1296 for (
unsigned int edge = 0; edge < n_edges; ++edge)
1299 in >> dummy >> dummy;
1305 skip_empty_lines(in);
1314 AssertThrow(line ==
"CrackedEdges", ExcInvalidDBMESHInput(line));
1317 for (
unsigned int edge = 0; edge < n_edges; ++edge)
1320 in >> dummy >> dummy;
1326 skip_empty_lines(in);
1332 AssertThrow(line ==
"Quadrilaterals", ExcInvalidDBMESHInput(line));
1334 static constexpr std::array<unsigned int, 8> local_vertex_numbering = {
1335 {0, 1, 5, 4, 2, 3, 7, 6}};
1336 std::vector<CellData<dim>> cells;
1338 unsigned int n_cells;
1340 for (
unsigned int cell = 0; cell < n_cells; ++cell)
1344 cells.emplace_back();
1348 cells.back().vertices[dim == 3 ? local_vertex_numbering[i] :
1352 (
static_cast<unsigned int>(cells.back().vertices[i]) <=
1354 ExcInvalidVertexIndex(cell, cells.back().vertices[i]));
1356 --cells.back().vertices[i];
1364 skip_empty_lines(in);
1372 while (getline(in, line), ((line.find(
"End") == std::string::npos) && (in)))
1378 apply_grid_fixup_functions(
vertices, cells, subcelldata);
1384template <
int dim,
int spacedim>
1388 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
1391 const auto reference_cell = ReferenceCells::get_hypercube<dim>();
1395 std::getline(in, line);
1397 unsigned int n_vertices;
1398 unsigned int n_cells;
1402 std::getline(in, line);
1405 std::getline(in, line);
1408 for (
unsigned int i = 0; i < 8; ++i)
1409 std::getline(in, line);
1412 std::vector<CellData<dim>> cells(n_cells);
1423 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
1424 in >> cell.vertices[reference_cell.exodusii_vertex_to_deal_vertex(i)];
1428 std::vector<Point<spacedim>>
vertices(n_vertices);
1431 for (
unsigned int d = 0; d < spacedim; ++d)
1433 for (
unsigned int d = spacedim; d < 3; ++d)
1442 apply_grid_fixup_functions(
vertices, cells, subcelldata);
1448template <
int dim,
int spacedim>
1452 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
1465 std::stringstream whole_file;
1470 std::getline(in, line);
1478 if ((line.size() > 0) && (line.back() ==
'\r'))
1479 line.erase(line.size() - 1);
1484 if (line.find(
'#') != std::string::npos)
1485 line.erase(line.find(
'#'), std::string::npos);
1486 while ((line.size() > 0) && (line.back() ==
' '))
1487 line.erase(line.size() - 1);
1489 if (line.size() > 0)
1490 whole_file <<
'\n' << line;
1509 unsigned int version_major, version_minor;
1510 whole_file >> version_major >> version_minor;
1511 AssertThrow((version_major == 0) && (version_minor == 1),
1512 ExcMessage(
"deal.II can currently only read version 0.1 "
1513 "of the mphtxt file format."));
1518 unsigned int n_tags;
1519 whole_file >> n_tags;
1520 for (
unsigned int i = 0; i < n_tags; ++i)
1523 while (whole_file.peek() ==
'\n')
1525 std::getline(whole_file, dummy);
1531 unsigned int n_types;
1532 whole_file >> n_types;
1533 for (
unsigned int i = 0; i < n_types; ++i)
1536 while (whole_file.peek() ==
'\n')
1538 std::getline(whole_file, dummy);
1558 whole_file >> dummy >> dummy >> dummy;
1562 while (whole_file.peek() ==
'\n')
1564 std::getline(whole_file, s);
1566 ExcMessage(
"Expected '4 Mesh', but got '" + s +
"'."));
1569 unsigned int version;
1570 whole_file >> version;
1574 unsigned int file_space_dim;
1575 whole_file >> file_space_dim;
1579 "The mesh file uses a different number of space dimensions "
1580 "than the triangulation you want to read it into."));
1582 unsigned int n_vertices;
1583 whole_file >> n_vertices;
1585 unsigned int starting_vertex_index;
1586 whole_file >> starting_vertex_index;
1588 std::vector<Point<spacedim>>
vertices(n_vertices);
1589 for (
unsigned int v = 0; v < n_vertices; ++v)
1620 std::vector<CellData<dim>> cells;
1623 unsigned int n_types;
1624 whole_file >> n_types;
1625 for (
unsigned int type = 0; type < n_types; ++type)
1632 whole_file >> dummy;
1636 std::string object_name;
1637 whole_file >> object_name;
1639 static const std::map<std::string, ReferenceCell> name_to_type = {
1649 AssertThrow(name_to_type.find(object_name) != name_to_type.end(),
1650 ExcMessage(
"The input file contains a cell type <" +
1652 "> that the reader does not "
1653 "current support"));
1654 const ReferenceCell object_type = name_to_type.at(object_name);
1656 unsigned int n_vertices_per_element;
1657 whole_file >> n_vertices_per_element;
1659 unsigned int n_elements;
1660 whole_file >> n_elements;
1674 ExcMessage(
"Triangles should not appear in input files "
1682 "Quadrilaterals should not appear in input files "
1689 ExcMessage(
"Tetrahedra should not appear in input files "
1690 "for 1d or 2d meshes."));
1697 "Prisms (wedges) should not appear in input files "
1698 "for 1d or 2d meshes."));
1717 std::vector<unsigned int> vertices_for_this_element(
1718 n_vertices_per_element);
1719 for (
unsigned int e = 0; e < n_elements; ++e)
1722 for (
unsigned int v = 0; v < n_vertices_per_element; ++v)
1724 whole_file >> vertices_for_this_element[v];
1725 vertices_for_this_element[v] -= starting_vertex_index;
1734 cells.emplace_back();
1735 cells.back().vertices = vertices_for_this_element;
1741 vertices_for_this_element;
1749 cells.emplace_back();
1750 cells.back().vertices = vertices_for_this_element;
1756 vertices_for_this_element;
1764 cells.emplace_back();
1765 cells.back().vertices = vertices_for_this_element;
1778 unsigned int n_geom_entity_indices;
1779 whole_file >> n_geom_entity_indices;
1781 (n_geom_entity_indices == n_elements),
1788 if (n_geom_entity_indices != 0)
1790 for (
unsigned int e = 0; e < n_geom_entity_indices; ++e)
1793 unsigned int geometric_entity_index;
1794 whole_file >> geometric_entity_index;
1800 cells[cells.size() - n_elements + e].material_id =
1801 geometric_entity_index;
1806 .boundary_id = geometric_entity_index;
1812 cells[cells.size() - n_elements + e].material_id =
1813 geometric_entity_index;
1818 .boundary_id = geometric_entity_index;
1824 cells[cells.size() - n_elements + e].material_id =
1825 geometric_entity_index;
1851 if (line.vertices[1] < line.vertices[0])
1852 std::swap(line.vertices[0], line.vertices[1]);
1857 return std::lexicographical_compare(a.vertices.begin(),
1879 Assert((face.vertices.size() == 3) || (face.vertices.size() == 4),
1881 std::sort(face.vertices.begin(), face.vertices.end());
1886 return std::lexicographical_compare(a.vertices.begin(),
1897 for (
const auto &face : cell->face_iterators())
1898 if (face->at_boundary())
1904 std::array<unsigned int, 2> face_vertex_indices = {
1905 {face->vertex_index(0), face->vertex_index(1)}};
1906 if (face_vertex_indices[0] > face_vertex_indices[1])
1907 std::swap(face_vertex_indices[0], face_vertex_indices[1]);
1913 face_vertex_indices,
1915 const std::array<unsigned int, 2>
1916 &face_vertex_indices) ->
bool {
1917 return std::lexicographical_compare(
1920 face_vertex_indices.begin(),
1921 face_vertex_indices.end());
1925 (p->vertices[0] == face_vertex_indices[0]) &&
1926 (p->vertices[1] == face_vertex_indices[1]))
1928 face->set_boundary_id(p->boundary_id);
1936 std::vector<unsigned int> face_vertex_indices(
1937 face->n_vertices());
1938 for (
unsigned int v = 0; v < face->n_vertices(); ++v)
1939 face_vertex_indices[v] = face->vertex_index(v);
1940 std::sort(face_vertex_indices.begin(),
1941 face_vertex_indices.end());
1947 face_vertex_indices,
1949 const std::vector<unsigned int>
1950 &face_vertex_indices) ->
bool {
1951 return std::lexicographical_compare(
1954 face_vertex_indices.begin(),
1955 face_vertex_indices.end());
1959 (p->vertices == face_vertex_indices))
1961 face->set_boundary_id(p->boundary_id);
1966 for (
unsigned int e = 0; e < face->n_lines(); ++e)
1968 const auto edge = face->line(e);
1970 std::array<unsigned int, 2> edge_vertex_indices = {
1971 {edge->vertex_index(0), edge->vertex_index(1)}};
1972 if (edge_vertex_indices[0] > edge_vertex_indices[1])
1973 std::swap(edge_vertex_indices[0],
1974 edge_vertex_indices[1]);
1980 edge_vertex_indices,
1982 const std::array<unsigned int, 2>
1983 &edge_vertex_indices) ->
bool {
1984 return std::lexicographical_compare(
1987 edge_vertex_indices.begin(),
1988 edge_vertex_indices.end());
1992 (p->vertices[0] == edge_vertex_indices[0]) &&
1993 (p->vertices[1] == edge_vertex_indices[1]))
1995 edge->set_boundary_id(p->boundary_id);
2004template <
int dim,
int spacedim>
2008 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
2011 unsigned int n_vertices;
2012 unsigned int n_cells;
2018 std::array<std::map<int, int>, 4> tag_maps;
2023 unsigned int gmsh_file_format = 0;
2024 if (line ==
"@f$NOD")
2025 gmsh_file_format = 10;
2026 else if (line ==
"@f$MeshFormat")
2027 gmsh_file_format = 20;
2033 if (gmsh_file_format == 20)
2036 unsigned int file_type, data_size;
2038 in >> version >> file_type >> data_size;
2041 gmsh_file_format =
static_cast<unsigned int>(version * 10);
2049 AssertThrow(line ==
"@f$EndMeshFormat", ExcInvalidGMSHInput(line));
2053 if (line ==
"@f$PhysicalNames")
2059 while (line !=
"@f$EndPhysicalNames");
2064 if (line ==
"@f$Entities")
2066 unsigned long n_points, n_curves, n_surfaces, n_volumes;
2068 in >> n_points >> n_curves >> n_surfaces >> n_volumes;
2069 for (
unsigned int i = 0; i < n_points; ++i)
2073 unsigned int n_physicals;
2074 double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y,
2078 if (gmsh_file_format > 40)
2080 in >> tag >> box_min_x >> box_min_y >> box_min_z >>
2082 box_max_x = box_min_x;
2083 box_max_y = box_min_y;
2084 box_max_z = box_min_z;
2088 in >> tag >> box_min_x >> box_min_y >> box_min_z >>
2089 box_max_x >> box_max_y >> box_max_z >> n_physicals;
2094 ExcMessage(
"More than one tag is not supported!"));
2096 int physical_tag = 0;
2097 for (
unsigned int j = 0; j < n_physicals; ++j)
2099 tag_maps[0][tag] = physical_tag;
2101 for (
unsigned int i = 0; i < n_curves; ++i)
2105 unsigned int n_physicals;
2106 double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y,
2110 in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >>
2111 box_max_y >> box_max_z >> n_physicals;
2115 ExcMessage(
"More than one tag is not supported!"));
2117 int physical_tag = 0;
2118 for (
unsigned int j = 0; j < n_physicals; ++j)
2120 tag_maps[1][tag] = physical_tag;
2125 for (
unsigned int j = 0; j < n_points; ++j)
2129 for (
unsigned int i = 0; i < n_surfaces; ++i)
2133 unsigned int n_physicals;
2134 double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y,
2138 in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >>
2139 box_max_y >> box_max_z >> n_physicals;
2143 ExcMessage(
"More than one tag is not supported!"));
2145 int physical_tag = 0;
2146 for (
unsigned int j = 0; j < n_physicals; ++j)
2148 tag_maps[2][tag] = physical_tag;
2153 for (
unsigned int j = 0; j < n_curves; ++j)
2156 for (
unsigned int i = 0; i < n_volumes; ++i)
2160 unsigned int n_physicals;
2161 double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y,
2165 in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >>
2166 box_max_y >> box_max_z >> n_physicals;
2170 ExcMessage(
"More than one tag is not supported!"));
2172 int physical_tag = 0;
2173 for (
unsigned int j = 0; j < n_physicals; ++j)
2175 tag_maps[3][tag] = physical_tag;
2180 for (
unsigned int j = 0; j < n_surfaces; ++j)
2184 AssertThrow(line ==
"@f$EndEntities", ExcInvalidGMSHInput(line));
2189 if (line ==
"@f$PartitionedEntities")
2195 while (line !=
"@f$EndPartitionedEntities");
2202 AssertThrow(line ==
"@f$Nodes", ExcInvalidGMSHInput(line));
2206 int n_entity_blocks = 1;
2207 if (gmsh_file_format > 40)
2211 in >> n_entity_blocks >> n_vertices >> min_node_tag >> max_node_tag;
2213 else if (gmsh_file_format == 40)
2215 in >> n_entity_blocks >> n_vertices;
2219 std::vector<Point<spacedim>>
vertices(n_vertices);
2226 unsigned int global_vertex = 0;
2227 for (
int entity_block = 0; entity_block < n_entity_blocks; ++entity_block)
2230 unsigned long numNodes;
2232 if (gmsh_file_format < 40)
2234 numNodes = n_vertices;
2241 int tagEntity, dimEntity;
2242 in >> tagEntity >> dimEntity >> parametric >> numNodes;
2245 std::vector<int> vertex_numbers;
2247 if (gmsh_file_format > 40)
2248 for (
unsigned long vertex_per_entity = 0;
2249 vertex_per_entity < numNodes;
2250 ++vertex_per_entity)
2252 in >> vertex_number;
2253 vertex_numbers.push_back(vertex_number);
2256 for (
unsigned long vertex_per_entity = 0; vertex_per_entity < numNodes;
2257 ++vertex_per_entity, ++global_vertex)
2263 if (gmsh_file_format > 40)
2265 vertex_number = vertex_numbers[vertex_per_entity];
2266 in >> x[0] >> x[1] >> x[2];
2269 in >> vertex_number >> x[0] >> x[1] >> x[2];
2271 for (
unsigned int d = 0; d < spacedim; ++d)
2277 if (parametric != 0)
2292 static const std::string end_nodes_marker[] = {
"@f$ENDNOD",
"@f$EndNodes"};
2293 AssertThrow(line == end_nodes_marker[gmsh_file_format == 10 ? 0 : 1],
2294 ExcInvalidGMSHInput(line));
2298 static const std::string begin_elements_marker[] = {
"@f$ELM",
"@f$Elements"};
2299 AssertThrow(line == begin_elements_marker[gmsh_file_format == 10 ? 0 : 1],
2300 ExcInvalidGMSHInput(line));
2303 if (gmsh_file_format > 40)
2307 in >> n_entity_blocks >> n_cells >> min_node_tag >> max_node_tag;
2309 else if (gmsh_file_format == 40)
2311 in >> n_entity_blocks >> n_cells;
2315 n_entity_blocks = 1;
2322 std::vector<CellData<dim>> cells;
2324 std::map<unsigned int, types::boundary_id> boundary_ids_1d;
2330 std::map<unsigned int, unsigned int> vertex_counts;
2333 static constexpr std::array<unsigned int, 8> local_vertex_numbering = {
2334 {0, 1, 5, 4, 2, 3, 7, 6}};
2335 unsigned int global_cell = 0;
2336 for (
int entity_block = 0; entity_block < n_entity_blocks; ++entity_block)
2338 unsigned int material_id;
2339 unsigned long numElements;
2342 if (gmsh_file_format < 40)
2346 numElements = n_cells;
2348 else if (gmsh_file_format == 40)
2350 int tagEntity, dimEntity;
2351 in >> tagEntity >> dimEntity >> cell_type >> numElements;
2352 material_id = tag_maps[dimEntity][tagEntity];
2357 int tagEntity, dimEntity;
2358 in >> dimEntity >> tagEntity >> cell_type >> numElements;
2359 material_id = tag_maps[dimEntity][tagEntity];
2362 for (
unsigned int cell_per_entity = 0; cell_per_entity < numElements;
2363 ++cell_per_entity, ++global_cell)
2372 unsigned int nod_num;
2393 unsigned int elm_number = 0;
2394 if (gmsh_file_format < 40)
2400 if (gmsh_file_format < 20)
2406 else if (gmsh_file_format < 40)
2411 unsigned int n_tags;
2418 for (
unsigned int i = 1; i < n_tags; ++i)
2423 else if (cell_type == 2)
2425 else if (cell_type == 3)
2427 else if (cell_type == 4)
2429 else if (cell_type == 5)
2440 else if (cell_type == 2)
2442 else if (cell_type == 3)
2444 else if (cell_type == 4)
2446 else if (cell_type == 5)
2472 if (((cell_type == 1) && (dim == 1)) ||
2473 ((cell_type == 2) && (dim == 2)) ||
2474 ((cell_type == 3) && (dim == 2)) ||
2475 ((cell_type == 4) && (dim == 3)) ||
2476 ((cell_type == 5) && (dim == 3)))
2479 unsigned int vertices_per_cell = 0;
2481 vertices_per_cell = 2;
2482 else if (cell_type == 2)
2483 vertices_per_cell = 3;
2484 else if (cell_type == 3)
2485 vertices_per_cell = 4;
2486 else if (cell_type == 4)
2487 vertices_per_cell = 4;
2488 else if (cell_type == 5)
2489 vertices_per_cell = 8;
2493 "Number of nodes does not coincide with the "
2494 "number required for this object"));
2497 cells.emplace_back();
2499 cell.
vertices.resize(vertices_per_cell);
2500 for (
unsigned int i = 0; i < vertices_per_cell; ++i)
2503 if (vertices_per_cell ==
2506 local_vertex_numbering[i] :
2514 std::numeric_limits<types::material_id>::max(),
2518 std::numeric_limits<types::material_id>::max()));
2526 for (
unsigned int i = 0; i < vertices_per_cell; ++i)
2530 ExcInvalidVertexIndexGmsh(cell_per_entity,
2536 vertex_counts[vertex] += 1u;
2540 else if ((cell_type == 1) &&
2541 ((dim == 2) || (dim == 3)))
2550 std::numeric_limits<types::boundary_id>::max(),
2554 std::numeric_limits<types::boundary_id>::max()));
2565 for (
unsigned int &vertex :
2574 ExcInvalidVertexIndex(cell_per_entity,
2579 else if ((cell_type == 2 || cell_type == 3) &&
2583 unsigned int vertices_per_cell = 0;
2586 vertices_per_cell = 3;
2587 else if (cell_type == 3)
2588 vertices_per_cell = 4;
2596 for (
unsigned int i = 0; i < vertices_per_cell; ++i)
2601 std::numeric_limits<types::boundary_id>::max(),
2605 std::numeric_limits<types::boundary_id>::max()));
2616 for (
unsigned int &vertex :
2625 ExcInvalidVertexIndex(cell_per_entity, vertex));
2629 else if (cell_type == 15)
2632 unsigned int node_index = 0;
2633 if (gmsh_file_format < 20)
2653 AssertThrow(
false, ExcGmshUnsupportedGeometry(cell_type));
2661 static const std::string end_elements_marker[] = {
"@f$ENDELM",
"@f$EndElements"};
2662 AssertThrow(line == end_elements_marker[gmsh_file_format == 10 ? 0 : 1],
2663 ExcInvalidGMSHInput(line));
2676 for (
const auto &pair : vertex_counts)
2677 if (pair.second == 1u)
2678 boundary_id_pairs.emplace_back(
vertices[pair.first],
2679 boundary_ids_1d[pair.first]);
2681 apply_grid_fixup_functions(
vertices, cells, subcelldata);
2687 assign_1d_boundary_ids(boundary_id_pairs, *tria);
2692#ifdef DEAL_II_GMSH_WITH_API
2693template <
int dim,
int spacedim>
2697 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
2699 const std::map<int, std::uint8_t> gmsh_to_dealii_type = {
2700 {{15, 0}, {1, 1}, {2, 2}, {3, 3}, {4, 4}, {7, 5}, {6, 6}, {5, 7}}};
2703 const std::array<std::vector<unsigned int>, 8> gmsh_to_dealii = {
2710 {{0, 1, 2, 3, 4, 5}},
2711 {{0, 1, 3, 2, 4, 5, 7, 6}}}};
2713 std::vector<Point<spacedim>>
vertices;
2714 std::vector<CellData<dim>> cells;
2716 std::map<unsigned int, types::boundary_id> boundary_ids_1d;
2722 std::map<unsigned int, unsigned int> vertex_counts;
2725 gmsh::option::setNumber(
"General.Verbosity", 0);
2729 ExcMessage(
"You are trying to read a gmsh file with dimension " +
2730 std::to_string(gmsh::model::getDimension()) +
2731 " into a grid of dimension " + std::to_string(dim)));
2736 gmsh::model::mesh::removeDuplicateNodes();
2737 gmsh::model::mesh::renumberNodes();
2738 std::vector<std::size_t> node_tags;
2739 std::vector<double> coord;
2740 std::vector<double> parametricCoord;
2741 gmsh::model::mesh::getNodes(
2742 node_tags, coord, parametricCoord, -1, -1,
false,
false);
2744 for (
unsigned int i = 0; i < node_tags.size(); ++i)
2748 for (
unsigned int d = 0; d < spacedim; ++d)
2752 for (
unsigned int d = spacedim; d < 3; ++d)
2754 ExcMessage(
"The grid you are reading contains nodes that are "
2755 "nonzero in the coordinate with index " +
2757 ", but you are trying to save "
2758 "it on a grid embedded in a " +
2759 std::to_string(spacedim) +
" dimensional space."));
2766 std::vector<std::pair<int, int>> entities;
2767 gmsh::model::getEntities(entities);
2769 for (
const auto &e : entities)
2772 const int &entity_dim = e.first;
2773 const int &entity_tag = e.second;
2779 std::vector<int> physical_tags;
2780 gmsh::model::getPhysicalGroupsForEntity(entity_dim,
2785 if (physical_tags.size())
2786 for (
auto physical_tag : physical_tags)
2789 gmsh::model::getPhysicalName(entity_dim, physical_tag, name);
2796 std::map<std::string, int> id_names;
2798 bool found_unrecognized_tag =
false;
2799 bool found_boundary_id =
false;
2802 for (
const auto &it : id_names)
2804 const auto &name = it.first;
2805 const auto &
id = it.second;
2806 if (entity_dim == dim && name ==
"MaterialID")
2809 found_boundary_id =
true;
2811 else if (entity_dim < dim && name ==
"BoundaryID")
2814 found_boundary_id =
true;
2816 else if (name ==
"ManifoldID")
2822 found_unrecognized_tag =
true;
2827 if (found_unrecognized_tag && !found_boundary_id)
2828 boundary_id = physical_tag;
2836 boundary_id = physical_tag;
2842 std::vector<int> element_types;
2843 std::vector<std::vector<std::size_t>> element_ids, element_nodes;
2844 gmsh::model::mesh::getElements(
2845 element_types, element_ids, element_nodes, entity_dim, entity_tag);
2847 for (
unsigned int i = 0; i < element_types.size(); ++i)
2849 const auto &type = gmsh_to_dealii_type.at(element_types[i]);
2850 const auto n_vertices = gmsh_to_dealii[type].size();
2851 const auto &elements = element_ids[i];
2852 const auto &nodes = element_nodes[i];
2853 for (
unsigned int j = 0; j < elements.size(); ++j)
2855 if (entity_dim == dim)
2857 cells.emplace_back(n_vertices);
2858 auto &cell = cells.back();
2859 for (
unsigned int v = 0; v < n_vertices; ++v)
2862 nodes[n_vertices * j + gmsh_to_dealii[type][v]] - 1;
2864 vertex_counts[cell.vertices[v]] += 1u;
2866 cell.manifold_id = manifold_id;
2867 cell.material_id = boundary_id;
2869 else if (entity_dim == 2)
2873 for (
unsigned int v = 0; v < n_vertices; ++v)
2875 nodes[n_vertices * j + gmsh_to_dealii[type][v]] - 1;
2877 face.manifold_id = manifold_id;
2878 face.boundary_id = boundary_id;
2880 else if (entity_dim == 1)
2884 for (
unsigned int v = 0; v < n_vertices; ++v)
2886 nodes[n_vertices * j + gmsh_to_dealii[type][v]] - 1;
2888 line.manifold_id = manifold_id;
2889 line.boundary_id = boundary_id;
2891 else if (entity_dim == 0)
2895 for (
unsigned int j = 0; j < elements.size(); ++j)
2896 boundary_ids_1d[nodes[j] - 1] = boundary_id;
2907 for (
const auto &pair : vertex_counts)
2908 if (pair.second == 1u)
2909 boundary_id_pairs.emplace_back(
vertices[pair.first],
2910 boundary_ids_1d[pair.first]);
2912 apply_grid_fixup_functions(
vertices, cells, subcelldata);
2918 assign_1d_boundary_ids(boundary_id_pairs, *tria);
2927template <
int dim,
int spacedim>
2930 std::string &header,
2931 std::vector<unsigned int> &tecplot2deal,
2932 unsigned int &n_vars,
2933 unsigned int &n_vertices,
2934 unsigned int &n_cells,
2935 std::vector<unsigned int> &IJK,
2960 std::transform(header.begin(), header.end(), header.begin(), ::toupper);
2964 std::replace(header.begin(), header.end(),
'\t',
' ');
2965 std::replace(header.begin(), header.end(),
',',
' ');
2966 std::replace(header.begin(), header.end(),
'\n',
' ');
2970 std::string::size_type pos = header.find(
'=');
2972 while (pos !=
static_cast<std::string::size_type
>(std::string::npos))
2973 if (header[pos + 1] ==
' ')
2974 header.erase(pos + 1, 1);
2975 else if (header[pos - 1] ==
' ')
2977 header.erase(pos - 1, 1);
2981 pos = header.find(
'=', ++pos);
2984 std::vector<std::string> entries =
2988 for (
unsigned int i = 0; i < entries.size(); ++i)
2997 tecplot2deal[0] = 0;
3000 while (entries[i][0] ==
'"')
3002 if (entries[i] ==
"\"X\"")
3003 tecplot2deal[0] = n_vars;
3004 else if (entries[i] ==
"\"Y\"")
3010 tecplot2deal[1] = n_vars;
3012 else if (entries[i] ==
"\"Z\"")
3018 tecplot2deal[2] = n_vars;
3030 "Tecplot file must contain at least one variable for each dimension"));
3031 for (
unsigned int d = 1; d < dim; ++d)
3033 tecplot2deal[d] > 0,
3035 "Tecplot file must contain at least one variable for each dimension."));
3040 "ZONETYPE=FELINESEG") &&
3044 "ZONETYPE=FEQUADRILATERAL") &&
3048 "ZONETYPE=FEBRICK") &&
3056 "The tecplot file contains an unsupported ZONETYPE."));
3059 "DATAPACKING=POINT"))
3062 "DATAPACKING=BLOCK"))
3085 "ET=QUADRILATERAL") &&
3097 "The tecplot file contains an unsupported ElementType."));
3105 dim > 1 || IJK[1] == 1,
3107 "Parameter 'J=' found in tecplot, although this is only possible for dimensions greater than 1."));
3113 dim > 2 || IJK[2] == 1,
3115 "Parameter 'K=' found in tecplot, although this is only possible for dimensions greater than 2."));
3130 for (
unsigned int d = 0; d < dim; ++d)
3135 "Tecplot file does not contain a complete and consistent set of parameters"));
3136 n_vertices *= IJK[d];
3137 n_cells *= (IJK[d] - 1);
3145 "Tecplot file does not contain a complete and consistent set of parameters"));
3151 n_cells = *std::max_element(IJK.begin(), IJK.end());
3155 "Tecplot file does not contain a complete and consistent set of parameters"));
3165 const unsigned int dim = 2;
3166 const unsigned int spacedim = 2;
3167 Assert(tria !=
nullptr, ExcNoTriangulationSelected());
3171 skip_comment_lines(in,
'#');
3174 std::string line, header;
3181 std::string letters =
"abcdfghijklmnopqrstuvwxyzABCDFGHIJKLMNOPQRSTUVWXYZ";
3184 while (line.find_first_of(letters) != std::string::npos)
3186 header +=
" " + line;
3193 std::vector<unsigned int> tecplot2deal(dim);
3194 std::vector<unsigned int> IJK(dim);
3195 unsigned int n_vars, n_vertices, n_cells;
3196 bool structured, blocked;
3198 parse_tecplot_header(header,
3213 std::vector<Point<spacedim>>
vertices(n_vertices + 1);
3216 std::vector<CellData<dim>> cells(n_cells);
3232 unsigned int next_index = 0;
3236 if (tecplot2deal[0] == 0)
3240 std::vector<std::string> first_var =
3243 for (
unsigned int i = 1; i < first_var.size() + 1; ++i)
3244 vertices[i][0] = std::strtod(first_var[i - 1].c_str(), &endptr);
3249 for (
unsigned int j = first_var.size() + 1; j < n_vertices + 1; ++j)
3257 for (
unsigned int i = 1; i < n_vars; ++i)
3266 if (next_index == dim && structured)
3269 if ((next_index < dim) && (i == tecplot2deal[next_index]))
3272 for (
unsigned int j = 1; j < n_vertices + 1; ++j)
3280 for (
unsigned int j = 1; j < n_vertices + 1; ++j)
3292 std::vector<double> vars(n_vars);
3297 std::vector<std::string> first_vertex =
3300 for (
unsigned int d = 0; d < dim; ++d)
3302 std::strtod(first_vertex[tecplot2deal[d]].c_str(), &endptr);
3306 for (
unsigned int v = 2; v < n_vertices + 1; ++v)
3308 for (
unsigned int i = 0; i < n_vars; ++i)
3314 for (
unsigned int i = 0; i < dim; ++i)
3315 vertices[v][i] = vars[tecplot2deal[i]];
3323 unsigned int I = IJK[0], J = IJK[1];
3325 unsigned int cell = 0;
3327 for (
unsigned int j = 0; j < J - 1; ++j)
3328 for (
unsigned int i = 1; i < I; ++i)
3330 cells[cell].vertices[0] = i + j * I;
3331 cells[cell].vertices[1] = i + 1 + j * I;
3332 cells[cell].vertices[2] = i + (j + 1) * I;
3333 cells[cell].vertices[3] = i + 1 + (j + 1) * I;
3337 std::vector<unsigned int> boundary_vertices(2 * I + 2 * J - 4);
3339 for (
unsigned int i = 1; i < I + 1; ++i)
3341 boundary_vertices[k] = i;
3343 boundary_vertices[k] = i + (J - 1) * I;
3346 for (
unsigned int j = 1; j < J - 1; ++j)
3348 boundary_vertices[k] = 1 + j * I;
3350 boundary_vertices[k] = I + j * I;
3369 for (
unsigned int i = 0; i < n_cells; ++i)
3386 apply_grid_fixup_functions(
vertices, cells, subcelldata);
3392template <
int dim,
int spacedim>
3401template <
int dim,
int spacedim>
3404 const unsigned int mesh_index,
3405 const bool remove_duplicates,
3407 const bool ignore_unsupported_types)
3409#ifdef DEAL_II_WITH_ASSIMP
3414 Assimp::Importer importer;
3417 const aiScene *scene =
3418 importer.ReadFile(filename.c_str(),
3419 aiProcess_RemoveComponent |
3420 aiProcess_JoinIdenticalVertices |
3421 aiProcess_ImproveCacheLocality | aiProcess_SortByPType |
3422 aiProcess_OptimizeGraph | aiProcess_OptimizeMeshes);
3428 ExcMessage(
"Input file contains no meshes."));
3431 (mesh_index < scene->mNumMeshes),
3434 unsigned int start_mesh =
3436 unsigned int end_mesh =
3441 std::vector<Point<spacedim>>
vertices;
3442 std::vector<CellData<dim>> cells;
3446 unsigned int v_offset = 0;
3447 unsigned int c_offset = 0;
3449 static constexpr std::array<unsigned int, 8> local_vertex_numbering = {
3450 {0, 1, 5, 4, 2, 3, 7, 6}};
3452 for (
unsigned int m = start_mesh; m < end_mesh; ++m)
3454 const aiMesh *mesh = scene->mMeshes[m];
3458 if ((dim == 2) && mesh->mPrimitiveTypes != aiPrimitiveType_POLYGON)
3461 ExcMessage(
"Incompatible mesh " + std::to_string(m) +
3462 "/" + std::to_string(scene->mNumMeshes)));
3465 else if ((dim == 1) && mesh->mPrimitiveTypes != aiPrimitiveType_LINE)
3468 ExcMessage(
"Incompatible mesh " + std::to_string(m) +
3469 "/" + std::to_string(scene->mNumMeshes)));
3473 const unsigned int n_vertices = mesh->mNumVertices;
3474 const aiVector3D *mVertices = mesh->mVertices;
3477 const unsigned int n_faces = mesh->mNumFaces;
3478 const aiFace *mFaces = mesh->mFaces;
3480 vertices.resize(v_offset + n_vertices);
3481 cells.resize(c_offset + n_faces);
3483 for (
unsigned int i = 0; i < n_vertices; ++i)
3484 for (
unsigned int d = 0; d < spacedim; ++d)
3485 vertices[i + v_offset][d] = mVertices[i][d];
3487 unsigned int valid_cell = c_offset;
3488 for (
unsigned int i = 0; i < n_faces; ++i)
3495 .vertices[dim == 3 ? local_vertex_numbering[f] :
3497 mFaces[i].mIndices[f] + v_offset;
3499 cells[valid_cell].material_id = m;
3505 ExcMessage(
"Face " + std::to_string(i) +
" of mesh " +
3506 std::to_string(m) +
" has " +
3507 std::to_string(mFaces[i].mNumIndices) +
3508 " vertices. We expected only " +
3513 cells.resize(valid_cell);
3518 v_offset += n_vertices;
3519 c_offset = valid_cell;
3526 if (remove_duplicates)
3531 unsigned int n_verts = 0;
3535 std::vector<unsigned int> considered_vertices;
3537 vertices, cells, subcelldata, considered_vertices, tol);
3541 apply_grid_fixup_functions(
vertices, cells, subcelldata);
3547 (void)remove_duplicates;
3549 (void)ignore_unsupported_types;
3554#ifdef DEAL_II_TRILINOS_WITH_SEACAS
3561 exodusii_name_to_type(
const std::string &type_name,
3562 const int n_nodes_per_element)
3568 std::string type_name_2 = type_name;
3569 std::transform(type_name_2.begin(),
3571 type_name_2.begin(),
3572 [](
unsigned char c) { return std::toupper(c); });
3573 const std::string
numbers =
"0123456789";
3574 type_name_2.erase(std::find_first_of(type_name_2.begin(),
3581 if (type_name_2 ==
"BAR" || type_name_2 ==
"BEAM" ||
3582 type_name_2 ==
"EDGE" || type_name_2 ==
"TRUSS")
3584 else if (type_name_2 ==
"TRI" || type_name_2 ==
"TRIANGLE")
3586 else if (type_name_2 ==
"QUAD" || type_name_2 ==
"QUADRILATERAL")
3588 else if (type_name_2 ==
"SHELL")
3590 if (n_nodes_per_element == 3)
3595 else if (type_name_2 ==
"TET" || type_name_2 ==
"TETRA" ||
3596 type_name_2 ==
"TETRAHEDRON")
3598 else if (type_name_2 ==
"PYRA" || type_name_2 ==
"PYRAMID")
3600 else if (type_name_2 ==
"WEDGE")
3602 else if (type_name_2 ==
"HEX" || type_name_2 ==
"HEXAHEDRON")
3612 template <
int dim,
int spacedim = dim>
3613 std::pair<SubCellData, std::vector<std::vector<int>>>
3614 read_exodusii_sidesets(
const int ex_id,
3615 const int n_side_sets,
3617 const bool apply_all_indicators_to_manifolds)
3620 std::vector<std::vector<int>> b_or_m_id_to_sideset_ids;
3622 b_or_m_id_to_sideset_ids.emplace_back();
3628 if (dim == spacedim && n_side_sets > 0)
3630 std::vector<int> side_set_ids(n_side_sets);
3631 int ierr = ex_get_ids(ex_id, EX_SIDE_SET, side_set_ids.data());
3632 AssertThrowExodusII(ierr);
3640 std::map<std::size_t, std::vector<int>> face_side_sets;
3641 for (
const int side_set_id : side_set_ids)
3644 int n_distribution_factors = -1;
3646 ierr = ex_get_set_param(ex_id,
3650 &n_distribution_factors);
3651 AssertThrowExodusII(ierr);
3654 std::vector<int> elements(n_sides);
3655 std::vector<int> faces(n_sides);
3656 ierr = ex_get_set(ex_id,
3661 AssertThrowExodusII(ierr);
3669 for (
int side_n = 0; side_n < n_sides; ++side_n)
3671 const long element_n = elements[side_n] - 1;
3672 const long face_n = faces[side_n] - 1;
3673 const std::size_t face_id =
3674 element_n * max_faces_per_cell + face_n;
3675 face_side_sets[face_id].push_back(side_set_id);
3681 std::vector<std::pair<std::size_t, std::vector<int>>>
3682 face_id_to_side_sets;
3683 face_id_to_side_sets.reserve(face_side_sets.size());
3684 for (
auto &pair : face_side_sets)
3687 face_id_to_side_sets.push_back(std::move(pair));
3691 std::sort(face_id_to_side_sets.begin(),
3692 face_id_to_side_sets.end(),
3693 [](
const auto &a,
const auto &b) {
3694 return std::lexicographical_compare(a.second.begin(),
3705 for (
const auto &pair : face_id_to_side_sets)
3707 const std::size_t face_id = pair.first;
3708 const std::vector<int> &face_sideset_ids = pair.second;
3709 if (face_sideset_ids != b_or_m_id_to_sideset_ids.back())
3714 ++current_b_or_m_id;
3715 b_or_m_id_to_sideset_ids.push_back(face_sideset_ids);
3716 Assert(current_b_or_m_id == b_or_m_id_to_sideset_ids.size() - 1,
3720 const unsigned int local_face_n = face_id % max_faces_per_cell;
3721 const CellData<dim> &cell = cells[face_id / max_faces_per_cell];
3724 const unsigned int deal_face_n =
3735 if (apply_all_indicators_to_manifolds)
3736 boundary_line.manifold_id = current_b_or_m_id;
3738 boundary_line.boundary_id = current_b_or_m_id;
3739 for (
unsigned int j = 0; j < face_reference_cell.
n_vertices();
3741 boundary_line.vertices[j] =
3743 deal_face_n, j, 0)];
3750 if (apply_all_indicators_to_manifolds)
3751 boundary_quad.manifold_id = current_b_or_m_id;
3753 boundary_quad.boundary_id = current_b_or_m_id;
3754 for (
unsigned int j = 0; j < face_reference_cell.
n_vertices();
3756 boundary_quad.vertices[j] =
3758 deal_face_n, j, 0)];
3765 return std::make_pair(std::move(subcelldata),
3766 std::move(b_or_m_id_to_sideset_ids));
3771template <
int dim,
int spacedim>
3774 const std::string &filename,
3775 const bool apply_all_indicators_to_manifolds)
3777#ifdef DEAL_II_TRILINOS_WITH_SEACAS
3779 int component_word_size =
sizeof(double);
3781 int floating_point_word_size = 0;
3782 float ex_version = 0.0;
3784 const int ex_id = ex_open(filename.c_str(),
3786 &component_word_size,
3787 &floating_point_word_size,
3790 ExcMessage(
"ExodusII failed to open the specified input file."));
3793 std::vector<char> string_temp(MAX_LINE_LENGTH + 1,
'\0');
3794 int mesh_dimension = 0;
3797 int n_element_blocks = 0;
3798 int n_node_sets = 0;
3799 int n_side_sets = 0;
3801 int ierr = ex_get_init(ex_id,
3809 AssertThrowExodusII(ierr);
3820 std::vector<Point<spacedim>>
vertices;
3823 std::vector<double> xs(n_nodes);
3824 std::vector<double> ys(n_nodes);
3825 std::vector<double> zs(n_nodes);
3827 ierr = ex_get_coord(ex_id, xs.data(), ys.data(), zs.data());
3828 AssertThrowExodusII(ierr);
3830 for (
int vertex_n = 0; vertex_n < n_nodes; ++vertex_n)
3835 vertices.emplace_back(xs[vertex_n]);
3838 vertices.emplace_back(xs[vertex_n], ys[vertex_n]);
3841 vertices.emplace_back(xs[vertex_n], ys[vertex_n], zs[vertex_n]);
3849 std::vector<int> element_block_ids(n_element_blocks);
3850 ierr = ex_get_ids(ex_id, EX_ELEM_BLOCK, element_block_ids.data());
3851 AssertThrowExodusII(ierr);
3853 std::vector<CellData<dim>> cells;
3854 cells.reserve(n_elements);
3858 for (
const int element_block_id : element_block_ids)
3860 std::fill(string_temp.begin(), string_temp.end(),
'\0');
3861 int n_block_elements = 0;
3862 int n_nodes_per_element = 0;
3863 int n_edges_per_element = 0;
3864 int n_faces_per_element = 0;
3865 int n_attributes_per_element = 0;
3868 ierr = ex_get_block(ex_id,
3873 &n_nodes_per_element,
3874 &n_edges_per_element,
3875 &n_faces_per_element,
3876 &n_attributes_per_element);
3877 AssertThrowExodusII(ierr);
3879 exodusii_name_to_type(string_temp.data(), n_nodes_per_element);
3882 "The ExodusII block " + std::to_string(element_block_id) +
3883 " with element type " + std::string(string_temp.data()) +
3885 ", which does not match the topological mesh dimension " +
3886 std::to_string(dim) +
"."));
3893 std::vector<int> connection(n_nodes_per_element * n_block_elements);
3894 ierr = ex_get_conn(ex_id,
3900 AssertThrowExodusII(ierr);
3902 for (
unsigned int elem_n = 0; elem_n < connection.size();
3903 elem_n += n_nodes_per_element)
3909 connection[elem_n + i] - 1;
3912 cells.push_back(std::move(cell));
3917 auto pair = read_exodusii_sidesets<dim, spacedim>(
3918 ex_id, n_side_sets, cells, apply_all_indicators_to_manifolds);
3919 ierr = ex_close(ex_id);
3920 AssertThrowExodusII(ierr);
3922 apply_grid_fixup_functions(
vertices, cells, pair.first);
3925 out.id_to_sideset_ids = std::move(pair.second);
3929 (void)apply_all_indicators_to_manifolds;
3936template <
int dim,
int spacedim>
3950 if (std::find_if(line.begin(), line.end(), [](
const char c) {
3955 for (
int i = line.size() - 1; i >= 0; --i)
3956 in.putback(line[i]);
3966template <
int dim,
int spacedim>
3969 const char comment_start)
3974 while (in.get(c) && c == comment_start)
3977 while (in.get() !=
'\n')
3987 skip_empty_lines(in);
3992template <
int dim,
int spacedim>
4010 double min_x =
vertices[cells[0].vertices[0]][0],
4011 max_x =
vertices[cells[0].vertices[0]][0],
4012 min_y =
vertices[cells[0].vertices[0]][1],
4013 max_y =
vertices[cells[0].vertices[0]][1];
4015 for (
unsigned int i = 0; i < cells.size(); ++i)
4017 for (
const auto vertex : cells[i].vertices)
4031 out <<
"# cell " << i << std::endl;
4033 for (
const auto vertex : cells[i].vertices)
4037 out <<
"set label \"" << i <<
"\" at " <<
center[0] <<
',' <<
center[1]
4038 <<
" center" << std::endl;
4041 for (
unsigned int f = 0; f < 2; ++f)
4045 <<
vertices[cells[i].vertices[(f + 1) % 4]][1] << std::endl;
4047 for (
unsigned int f = 2; f < 4; ++f)
4049 <<
',' <<
vertices[cells[i].vertices[(f + 1) % 4]][1] <<
" to "
4057 <<
"set nokey" << std::endl
4058 <<
"pl [" << min_x <<
':' << max_x <<
"][" << min_y <<
':' << max_y
4059 <<
"] " << min_y << std::endl
4060 <<
"pause -1" << std::endl;
4071 for (
const auto &cell : cells)
4138template <
int dim,
int spacedim>
4146 if (format == Default)
4147 name = search.
find(filename);
4149 name = search.
find(filename, default_suffix(format));
4152 if (format == Default)
4154 const std::string::size_type slashpos = name.find_last_of(
'/');
4155 const std::string::size_type dotpos = name.find_last_of(
'.');
4156 if (dotpos < name.size() &&
4157 (dotpos > slashpos || slashpos == std::string::npos))
4159 std::string ext = name.substr(dotpos + 1);
4160 format = parse_format(ext);
4164 if (format == assimp)
4168 else if (format == exodusii)
4170 read_exodusii(name);
4174 std::ifstream in(name);
4180template <
int dim,
int spacedim>
4184 if (format == Default)
4185 format = default_format;
4227 ExcMessage(
"There is no read_assimp(istream &) function. "
4228 "Use the read_assimp(string &filename, ...) "
4229 "functions, instead."));
4234 ExcMessage(
"There is no read_exodusii(istream &) function. "
4235 "Use the read_exodusii(string &filename, ...) "
4236 "function, instead."));
4247template <
int dim,
int spacedim>
4276 return ".unknown_format";
4282template <
int dim,
int spacedim>
4286 if (format_name ==
"dbmesh")
4289 if (format_name ==
"exodusii")
4292 if (format_name ==
"msh")
4295 if (format_name ==
"unv")
4298 if (format_name ==
"vtk")
4301 if (format_name ==
"vtu")
4305 if (format_name ==
"inp")
4308 if (format_name ==
"ucd")
4311 if (format_name ==
"xda")
4314 if (format_name ==
"tecplot")
4317 if (format_name ==
"dat")
4320 if (format_name ==
"plt")
4339template <
int dim,
int spacedim>
4343 return "dbmesh|exodusii|msh|unv|vtk|vtu|ucd|abaqus|xda|tecplot|assimp";
4350 template <
int dim,
int spacedim>
4351 Abaqus_to_UCD<dim, spacedim>::Abaqus_to_UCD()
4364 from_string(T &t,
const std::string &s, std::ios_base &(*f)(std::ios_base &))
4366 std::istringstream iss(s);
4367 return !(iss >> f >> t).fail();
4374 extract_int(
const std::string &s)
4377 for (
const char c : s)
4379 if (isdigit(c) != 0)
4386 from_string(number, tmp, std::dec);
4392 template <
int dim,
int spacedim>
4394 Abaqus_to_UCD<dim, spacedim>::read_in_abaqus(std::istream &input_stream)
4403 while (std::getline(input_stream, line))
4406 std::transform(line.begin(), line.end(), line.begin(), ::toupper);
4408 if (line.compare(
"*HEADING") == 0 || line.compare(0, 2,
"**") == 0 ||
4409 line.compare(0, 5,
"*PART") == 0)
4412 while (std::getline(input_stream, line))
4418 else if (line.compare(0, 5,
"*NODE") == 0)
4427 while (std::getline(input_stream, line))
4432 std::vector<double> node(spacedim + 1);
4434 std::istringstream iss(line);
4436 for (
unsigned int i = 0; i < spacedim + 1; ++i)
4437 iss >> node[i] >> comma;
4439 node_list.push_back(node);
4442 else if (line.compare(0, 8,
"*ELEMENT") == 0)
4457 const std::string before_material =
"ELSET=EB";
4458 const std::size_t idx = line.find(before_material);
4459 if (idx != std::string::npos)
4461 from_string(material,
4462 line.substr(idx + before_material.size()),
4468 while (std::getline(input_stream, line))
4473 std::istringstream iss(line);
4479 const unsigned int n_data_per_cell =
4481 std::vector<double> cell(n_data_per_cell);
4482 for (
unsigned int i = 0; i < n_data_per_cell; ++i)
4483 iss >> cell[i] >> comma;
4486 cell[0] =
static_cast<double>(material);
4487 cell_list.push_back(cell);
4490 else if (line.compare(0, 8,
"*SURFACE") == 0)
4501 const std::string name_key =
"NAME=";
4502 const std::size_t name_idx_start =
4503 line.find(name_key) + name_key.size();
4504 std::size_t name_idx_end = line.find(
',', name_idx_start);
4505 if (name_idx_end == std::string::npos)
4507 name_idx_end = line.size();
4509 const int b_indicator = extract_int(
4510 line.substr(name_idx_start, name_idx_end - name_idx_start));
4517 while (std::getline(input_stream, line))
4523 std::transform(line.begin(),
4531 std::istringstream iss(line);
4538 std::vector<double> quad_node_list;
4539 const std::string elset_name = line.substr(0, line.find(
','));
4540 if (elsets_list.count(elset_name) != 0)
4544 iss >> stmp >> temp >> face_number;
4546 const std::vector<int> cells = elsets_list[elset_name];
4547 for (
const int cell : cells)
4551 get_global_node_numbers(el_idx, face_number);
4552 quad_node_list.insert(quad_node_list.begin(),
4555 face_list.push_back(quad_node_list);
4562 iss >> el_idx >> comma >> temp >> face_number;
4564 get_global_node_numbers(el_idx, face_number);
4565 quad_node_list.insert(quad_node_list.begin(), b_indicator);
4567 face_list.push_back(quad_node_list);
4571 else if (line.compare(0, 6,
"*ELSET") == 0)
4575 std::string elset_name;
4577 const std::string elset_key =
"*ELSET, ELSET=";
4578 const std::size_t idx = line.find(elset_key);
4579 if (idx != std::string::npos)
4581 const std::string comma =
",";
4582 const std::size_t first_comma = line.find(comma);
4583 const std::size_t second_comma =
4584 line.find(comma, first_comma + 1);
4585 const std::size_t elset_name_start =
4586 line.find(elset_key) + elset_key.size();
4587 elset_name = line.substr(elset_name_start,
4588 second_comma - elset_name_start);
4598 std::vector<int> elements;
4599 const std::size_t generate_idx = line.find(
"GENERATE");
4600 if (generate_idx != std::string::npos)
4603 std::getline(input_stream, line);
4604 std::istringstream iss(line);
4614 iss >> elid_start >> comma >> elid_end;
4618 "While reading an ABAQUS file, the reader "
4619 "expected a comma but found a <") +
4620 comma +
"> in the line <" + line +
">."));
4622 elid_start <= elid_end,
4625 "While reading an ABAQUS file, the reader encountered "
4626 "a GENERATE statement in which the upper bound <") +
4628 "> for the element numbers is not larger or equal "
4629 "than the lower bound <" +
4633 if (iss.rdbuf()->in_avail() != 0)
4634 iss >> comma >> elis_step;
4638 "While reading an ABAQUS file, the reader "
4639 "expected a comma but found a <") +
4640 comma +
"> in the line <" + line +
">."));
4642 for (
int i = elid_start; i <= elid_end; i += elis_step)
4643 elements.push_back(i);
4644 elsets_list[elset_name] = elements;
4646 std::getline(input_stream, line);
4651 while (std::getline(input_stream, line))
4656 std::istringstream iss(line);
4661 iss >> elid >> comma;
4666 "While reading an ABAQUS file, the reader "
4667 "expected a comma but found a <") +
4668 comma +
"> in the line <" + line +
">."));
4670 elements.push_back(elid);
4674 elsets_list[elset_name] = elements;
4679 else if (line.compare(0, 5,
"*NSET") == 0)
4682 while (std::getline(input_stream, line))
4688 else if (line.compare(0, 14,
"*SOLID SECTION") == 0)
4692 const std::string elset_key =
"ELSET=";
4693 const std::size_t elset_start =
4694 line.find(
"ELSET=") + elset_key.size();
4695 const std::size_t elset_end = line.find(
',', elset_start + 1);
4696 const std::string elset_name =
4697 line.substr(elset_start, elset_end - elset_start);
4702 const std::string material_key =
"MATERIAL=";
4703 const std::size_t last_equal =
4704 line.find(
"MATERIAL=") + material_key.size();
4705 const std::size_t material_id_start = line.find(
'-', last_equal);
4707 from_string(material_id,
4708 line.substr(material_id_start + 1),
4712 const std::vector<int> &elset_cells = elsets_list[elset_name];
4713 for (
const int elset_cell : elset_cells)
4715 const int cell_id = elset_cell - 1;
4723 template <
int dim,
int spacedim>
4725 Abaqus_to_UCD<dim, spacedim>::get_global_node_numbers(
4726 const int face_cell_no,
4727 const int face_cell_face_no)
const
4737 if (face_cell_face_no == 1)
4739 quad_node_list[0] = cell_list[face_cell_no - 1][1];
4740 quad_node_list[1] = cell_list[face_cell_no - 1][2];
4742 else if (face_cell_face_no == 2)
4744 quad_node_list[0] = cell_list[face_cell_no - 1][2];
4745 quad_node_list[1] = cell_list[face_cell_no - 1][3];
4747 else if (face_cell_face_no == 3)
4749 quad_node_list[0] = cell_list[face_cell_no - 1][3];
4750 quad_node_list[1] = cell_list[face_cell_no - 1][4];
4752 else if (face_cell_face_no == 4)
4754 quad_node_list[0] = cell_list[face_cell_no - 1][4];
4755 quad_node_list[1] = cell_list[face_cell_no - 1][1];
4765 if (face_cell_face_no == 1)
4767 quad_node_list[0] = cell_list[face_cell_no - 1][1];
4768 quad_node_list[1] = cell_list[face_cell_no - 1][4];
4769 quad_node_list[2] = cell_list[face_cell_no - 1][3];
4770 quad_node_list[3] = cell_list[face_cell_no - 1][2];
4772 else if (face_cell_face_no == 2)
4774 quad_node_list[0] = cell_list[face_cell_no - 1][5];
4775 quad_node_list[1] = cell_list[face_cell_no - 1][8];
4776 quad_node_list[2] = cell_list[face_cell_no - 1][7];
4777 quad_node_list[3] = cell_list[face_cell_no - 1][6];
4779 else if (face_cell_face_no == 3)
4781 quad_node_list[0] = cell_list[face_cell_no - 1][1];
4782 quad_node_list[1] = cell_list[face_cell_no - 1][2];
4783 quad_node_list[2] = cell_list[face_cell_no - 1][6];
4784 quad_node_list[3] = cell_list[face_cell_no - 1][5];
4786 else if (face_cell_face_no == 4)
4788 quad_node_list[0] = cell_list[face_cell_no - 1][2];
4789 quad_node_list[1] = cell_list[face_cell_no - 1][3];
4790 quad_node_list[2] = cell_list[face_cell_no - 1][7];
4791 quad_node_list[3] = cell_list[face_cell_no - 1][6];
4793 else if (face_cell_face_no == 5)
4795 quad_node_list[0] = cell_list[face_cell_no - 1][3];
4796 quad_node_list[1] = cell_list[face_cell_no - 1][4];
4797 quad_node_list[2] = cell_list[face_cell_no - 1][8];
4798 quad_node_list[3] = cell_list[face_cell_no - 1][7];
4800 else if (face_cell_face_no == 6)
4802 quad_node_list[0] = cell_list[face_cell_no - 1][1];
4803 quad_node_list[1] = cell_list[face_cell_no - 1][5];
4804 quad_node_list[2] = cell_list[face_cell_no - 1][8];
4805 quad_node_list[3] = cell_list[face_cell_no - 1][4];
4818 return quad_node_list;
4821 template <
int dim,
int spacedim>
4823 Abaqus_to_UCD<dim, spacedim>::write_out_avs_ucd(std::ostream &output)
const
4832 const boost::io::ios_base_all_saver formatting_saver(output);
4836 output <<
"# Abaqus to UCD mesh conversion" << std::endl;
4837 output <<
"# Mesh type: AVS UCD" << std::endl;
4866 output << node_list.size() <<
"\t" << (cell_list.size() + face_list.size())
4867 <<
"\t0\t0\t0" << std::endl;
4870 output.precision(8);
4874 for (
const auto &node : node_list)
4877 output << node[0] <<
"\t";
4880 output.setf(std::ios::scientific, std::ios::floatfield);
4881 for (
unsigned int jj = 1; jj < spacedim + 1; ++jj)
4884 if (
std::abs(node[jj]) > tolerance)
4885 output << static_cast<double>(node[jj]) <<
"\t";
4887 output << 0.0 <<
"\t";
4890 output << 0.0 <<
"\t";
4892 output << std::endl;
4893 output.unsetf(std::ios::floatfield);
4897 for (
unsigned int ii = 0; ii < cell_list.size(); ++ii)
4899 output << ii + 1 <<
"\t" << cell_list[ii][0] <<
"\t"
4900 << (dim == 2 ?
"quad" :
"hex") <<
"\t";
4901 for (
unsigned int jj = 1; jj < GeometryInfo<dim>::vertices_per_cell + 1;
4903 output << cell_list[ii][jj] <<
"\t";
4905 output << std::endl;
4909 for (
unsigned int ii = 0; ii < face_list.size(); ++ii)
4911 output << ii + 1 <<
"\t" << face_list[ii][0] <<
"\t"
4912 << (dim == 2 ?
"line" :
"quad") <<
"\t";
4913 for (
unsigned int jj = 1; jj < GeometryInfo<dim>::vertices_per_face + 1;
4915 output << face_list[ii][jj] <<
"\t";
4917 output << std::endl;
4924#include "grid_in.inst"
void read_vtk(std::istream &in)
static void skip_empty_lines(std::istream &in)
void read_assimp(const std::string &filename, const unsigned int mesh_index=numbers::invalid_unsigned_int, const bool remove_duplicates=true, const double tol=1e-12, const bool ignore_unsupported_element_types=true)
void read_abaqus(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
static std::string default_suffix(const Format format)
void read_xda(std::istream &in)
void read_comsol_mphtxt(std::istream &in)
static void skip_comment_lines(std::istream &in, const char comment_start)
void attach_triangulation(Triangulation< dim, spacedim > &tria)
void read_msh(std::istream &in)
static Format parse_format(const std::string &format_name)
void read_vtu(std::istream &in)
void read_tecplot(std::istream &in)
ExodusIIData read_exodusii(const std::string &filename, const bool apply_all_indicators_to_manifolds=false)
void read_dbmesh(std::istream &in)
void read_ucd(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
void read(std::istream &in, Format format=Default)
static void debug_output_grid(const std::vector< CellData< dim > > &cells, const std::vector< Point< spacedim > > &vertices, std::ostream &out)
static std::string get_format_names()
void read_unv(std::istream &in)
static void parse_tecplot_header(std::string &header, std::vector< unsigned int > &tecplot2deal, unsigned int &n_vars, unsigned int &n_vertices, unsigned int &n_cells, std::vector< unsigned int > &IJK, bool &structured, bool &blocked)
std::string find(const std::string &filename, const char *open_mode="r")
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
unsigned int n_vertices() const
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const
unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
unsigned int get_dimension() const
ReferenceCell face_reference_cell(const unsigned int face_no) const
unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const
static ReferenceCell n_vertices_to_type(const int dim, const unsigned int n_vertices)
virtual void create_triangulation(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
void load(Archive &ar, const unsigned int version)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_FALLTHROUGH
IteratorRange< active_cell_iterator > active_cell_iterators() const
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcNeedsAssimp()
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
static ::ExceptionBase & ExcNeedsExodusII()
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcInvalidState()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
void consistently_order_cells(std::vector< CellData< dim > > &cells)
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
spacedim MeshType< dim - 1, spacedim > const std::set< types::boundary_id > & boundary_ids
constexpr const ReferenceCell Tetrahedron
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Pyramid
constexpr const ReferenceCell Invalid
constexpr const ReferenceCell Triangle
constexpr const ReferenceCell Hexahedron
constexpr const ReferenceCell Vertex
constexpr const ReferenceCell Line
std::pair< int, unsigned int > get_integer_at_position(const std::string &name, const unsigned int position)
std::vector< unsigned char > decode_base64(const std::string &base64_input)
std::vector< std::string > break_text_into_lines(const std::string &original_text, const unsigned int width, const char delimiter=' ')
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
bool match_at_string_start(const std::string &name, const std::string &pattern)
std::string decompress(const std::string &compressed_input)
const types::material_id invalid_material_id
const types::boundary_id internal_face_boundary_id
static const unsigned int invalid_unsigned_int
const types::manifold_id flat_manifold_id
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
std::vector< unsigned int > vertices
types::material_id material_id
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices()
std::vector< CellData< 2 > > boundary_quads
bool check_consistency(const unsigned int dim) const
std::vector< CellData< 1 > > boundary_lines