OGS
TetGenInterface.cpp
Go to the documentation of this file.
1
15#include "TetGenInterface.h"
16
17#include <cstddef>
18#include <fstream>
19#include <range/v3/numeric/accumulate.hpp>
20#include <string>
21
22#include "BaseLib/FileTools.h"
23#include "BaseLib/Logging.h"
24#include "BaseLib/StringTools.h"
25#include "GeoLib/Triangle.h"
28#include "MeshLib/Mesh.h"
29#include "MeshLib/Node.h"
31
32namespace FileIO
33{
35
36auto constructPointsFromNodes(std::vector<MeshLib::Node*> nodes)
37{
38 std::vector<GeoLib::Point*> points;
39 points.reserve(nodes.size());
40 std::transform(nodes.begin(), nodes.end(), std::back_inserter(points),
41 [](auto const* node)
42 { return new GeoLib::Point(*node, node->getID()); });
43 return points;
44}
45
46bool TetGenInterface::readTetGenGeometry(std::string const& geo_fname,
47 GeoLib::GEOObjects& geo_objects)
48{
49 std::ifstream poly_stream(geo_fname.c_str());
50
51 if (!poly_stream)
52 {
53 ERR("TetGenInterface::readTetGenGeometry() failed to open {:s}",
54 geo_fname);
55 return false;
56 }
57 std::string ext(BaseLib::getFileExtension(geo_fname));
58 if (ext != ".smesh")
59 {
60 ERR("TetGenInterface::readTetGenGeometry() - unknown file type (only "
61 "*.smesh is supported).");
62 return false;
63 }
64
65 std::vector<MeshLib::Node*> nodes;
66 if (!readNodesFromStream(poly_stream, nodes))
67 {
68 // remove nodes read until now
70 return false;
71 }
72 auto points = constructPointsFromNodes(nodes);
74
75 std::string geo_name(BaseLib::extractBaseNameWithoutExtension(geo_fname));
76 geo_objects.addPointVec(std::move(points), geo_name);
77 const std::vector<std::size_t>& id_map(
78 geo_objects.getPointVecObj(geo_name)->getIDMap());
79
80 std::vector<GeoLib::Surface*> surfaces;
81 if (!parseSmeshFacets(poly_stream, surfaces,
82 *geo_objects.getPointVec(geo_name), id_map))
83 {
84 // remove surfaces read until now but keep the points
86 }
87 geo_objects.addSurfaceVec(std::move(surfaces), geo_name,
89
90 return true;
91}
92
93std::size_t TetGenInterface::getNFacets(std::ifstream& input)
94{
95 std::string line;
96 while (!input.fail())
97 {
98 std::getline(input, line);
99 if (input.fail())
100 {
101 ERR("TetGenInterface::getNFacets(): Error reading number of "
102 "facets.");
103 return 0;
104 }
105
106 BaseLib::simplify(line);
107 if (line.empty() || line.compare(0, 1, "#") == 0)
108 {
109 continue;
110 }
111
112 const std::list<std::string> fields = BaseLib::splitString(line, ' ');
113 auto it = fields.begin();
114 const auto nFacets(BaseLib::str2number<std::size_t>(*it));
115 if (fields.size() > 1)
116 {
118 }
119 return nFacets;
120 }
121 return 0;
122}
123
125 std::ifstream& input,
126 std::vector<GeoLib::Surface*>& surfaces,
127 const std::vector<GeoLib::Point*>& points,
128 const std::vector<std::size_t>& pnt_id_map)
129{
130 const std::size_t nFacets(this->getNFacets(input));
131 std::string line;
132 surfaces.reserve(nFacets);
133 std::list<std::string>::const_iterator it;
134
135 const unsigned offset = (_zero_based_idx) ? 0 : 1;
136 std::vector<std::size_t> idx_map;
137
138 std::size_t k(0);
139 while (k < nFacets && !input.fail())
140 {
141 std::getline(input, line);
142 if (input.fail())
143 {
144 ERR("TetGenInterface::parseFacets(): Error reading facet {:d}.", k);
145 return false;
146 }
147
148 BaseLib::simplify(line);
149 if (line.empty() || line.compare(0, 1, "#") == 0)
150 {
151 continue;
152 }
153
154 // read facets
155 const std::list<std::string> point_fields =
156 BaseLib::splitString(line, ' ');
157 it = point_fields.begin();
158 const auto nPoints = BaseLib::str2number<std::size_t>(*it);
159 if (nPoints != 3)
160 {
161 ERR("Smesh-files are currently only supported for triangle "
162 "meshes.");
163 return false;
164 }
165 const std::size_t point_field_size =
166 (_boundary_markers) ? nPoints + 1 : nPoints;
167 if (point_fields.size() > point_field_size)
168 {
169 std::vector<std::size_t> point_ids;
170 for (std::size_t j(0); j < nPoints; ++j)
171 {
172 point_ids.push_back(
173 pnt_id_map[BaseLib::str2number<std::size_t>(*(++it)) -
174 offset]);
175 }
176
177 const std::size_t sfc_marker =
179 : 0;
180 const std::size_t idx =
181 std::find(idx_map.begin(), idx_map.end(), sfc_marker) -
182 idx_map.begin();
183 if (idx >= surfaces.size())
184 {
185 idx_map.push_back(sfc_marker);
186 surfaces.push_back(new GeoLib::Surface(points));
187 }
188 surfaces[idx]->addTriangle(point_ids[0], point_ids[1],
189 point_ids[2]);
190 }
191 else
192 {
193 ERR("TetGenInterface::parseFacets(): Error reading points for "
194 "facet {:d}.",
195 k);
196 return false;
197 }
198 ++k;
199 }
200 // here the poly-file potentially defines a number of points to mark holes
201 // within the volumes defined by the facets, these are ignored for now here
202 // the poly-file potentially defines a number of region attributes, these
203 // are ignored for now
204
205 auto const nTotalTriangles = ranges::accumulate(
206 surfaces, std::size_t{0}, {},
207 [](auto const* surface) { return surface->getNumberOfTriangles(); });
208 if (nTotalTriangles == nFacets)
209 {
210 return true;
211 }
212
213 ERR("TetGenInterface::parseFacets(): Number of expected total triangles "
214 "({:d}) does not match number of found triangles ({:d}).",
215 surfaces.size(),
216 nTotalTriangles);
217 return false;
218}
219
220MeshLib::Mesh* TetGenInterface::readTetGenMesh(std::string const& nodes_fname,
221 std::string const& ele_fname)
222{
223 std::ifstream ins_nodes(nodes_fname.c_str());
224 std::ifstream ins_ele(ele_fname.c_str());
225
226 if (!ins_nodes || !ins_ele)
227 {
228 if (!ins_nodes)
229 {
230 ERR("TetGenInterface::readTetGenMesh failed to open {:s}",
231 nodes_fname);
232 }
233 if (!ins_ele)
234 {
235 ERR("TetGenInterface::readTetGenMesh failed to open {:s}",
236 ele_fname);
237 }
238 return nullptr;
239 }
240
241 std::vector<MeshLib::Node*> nodes;
242 if (!readNodesFromStream(ins_nodes, nodes))
243 {
244 // remove nodes read until now
246 return nullptr;
247 }
248
249 std::vector<MeshLib::Element*> elements;
250 std::vector<int> materials;
251 if (!readElementsFromStream(ins_ele, elements, materials, nodes))
252 {
253 BaseLib::cleanupVectorElements(nodes, elements);
254 return nullptr;
255 }
256
257 MeshLib::Properties properties;
258 // Transmit material values if there is any material value != 0
259 if (std::any_of(materials.cbegin(), materials.cend(),
260 [](int m) { return m != 0; }))
261 {
262 auto* const mat_props = properties.createNewPropertyVector<int>(
263 "MaterialIDs", MeshLib::MeshItemType::Cell);
264 mat_props->assign(materials);
265 }
266
267 const std::string mesh_name(
269 return new MeshLib::Mesh(mesh_name, nodes, elements,
270 true /* compute_element_neighbors */, properties);
271}
272
274 std::vector<MeshLib::Node*>& nodes)
275{
276 std::string line;
277 std::getline(ins, line);
278 std::size_t n_nodes;
279 std::size_t dim;
280 std::size_t n_attributes;
281 bool boundary_markers;
282
283 while (!ins.fail())
284 {
285 BaseLib::simplify(line);
286 if (line.empty() || line.compare(0, 1, "#") == 0)
287 {
288 // this line is a comment - skip
289 std::getline(ins, line);
290 continue;
291 }
292 // read header line
293 bool header_okay = parseNodesFileHeader(line, n_nodes, dim,
294 n_attributes, boundary_markers);
295 if (!header_okay)
296 {
297 return false;
298 }
299 if (!parseNodes(ins, nodes, n_nodes, dim))
300 {
301 return false;
302 }
303 return true;
304 }
305 return false;
306}
307
308bool TetGenInterface::parseNodesFileHeader(std::string const& line,
309 std::size_t& n_nodes,
310 std::size_t& dim,
311 std::size_t& n_attributes,
312 bool& boundary_markers)
313{
314 std::list<std::string> pnt_header = BaseLib::splitString(line, ' ');
315 if (pnt_header.empty())
316 {
317 ERR("TetGenInterface::parseNodesFileHeader(): could not read number of "
318 "nodes specified in header.");
319 return false;
320 }
321 auto it = pnt_header.begin();
323 dim = (pnt_header.size() == 1) ? 3
325
326 if (pnt_header.size() < 4)
327 {
328 n_attributes = 0;
329 boundary_markers = false;
330 return true;
331 }
332
333 n_attributes = BaseLib::str2number<std::size_t>(*(++it));
334 boundary_markers = *(++it) == "1";
335
336 return true;
337}
338
339bool TetGenInterface::parseNodes(std::ifstream& ins,
340 std::vector<MeshLib::Node*>& nodes,
341 std::size_t n_nodes,
342 std::size_t dim)
343{
344 std::string line;
345 nodes.reserve(n_nodes);
346
347 std::size_t k(0);
348 while (k < n_nodes && !ins.fail())
349 {
350 std::vector<double> coordinates(dim);
351 std::getline(ins, line);
352 if (ins.fail())
353 {
354 ERR("TetGenInterface::parseNodes(): Error reading node {:d}.", k);
355 return false;
356 }
357
358 std::size_t id;
359 std::size_t pos_end = 0;
360 std::size_t pos_beg = line.find_first_not_of(' ', pos_end);
361 pos_end = line.find_first_of(" \n", pos_beg);
362
363 if (line.empty() || pos_beg == pos_end ||
364 line.compare(pos_beg, 1, "#") == 0)
365 {
366 continue;
367 }
368
369 if (pos_beg != std::string::npos && pos_end != std::string::npos)
370 {
372 line.substr(pos_beg, pos_end - pos_beg));
373 if (k == 0 && id == 0)
374 {
375 _zero_based_idx = true;
376 }
377 }
378 else
379 {
380 ERR("TetGenInterface::parseNodes(): Error reading ID of node {:d}.",
381 k);
382 return false;
383 }
384 // read coordinates
385 const unsigned offset = (_zero_based_idx) ? 0 : 1;
386 for (std::size_t i(0); i < dim; i++)
387 {
388 pos_beg = line.find_first_not_of(' ', pos_end);
389 pos_end = line.find_first_of(" \n", pos_beg);
390 if (pos_end == std::string::npos)
391 {
392 pos_end = line.size();
393 }
394 if (pos_beg != std::string::npos)
395 {
396 coordinates[i] = BaseLib::str2number<double>(
397 line.substr(pos_beg, pos_end - pos_beg));
398 }
399 else
400 {
401 ERR("TetGenInterface::parseNodes(): error reading coordinate "
402 "{:d} of node {:d}.",
403 i,
404 k);
405 return false;
406 }
407 }
408
409 nodes.push_back(new MeshLib::Node(coordinates.data(), id - offset));
410 // read attributes and boundary markers ... - at the moment we do not
411 // use this information
412 ++k;
413 }
414
415 return true;
416}
417
419 std::ifstream& ins,
420 std::vector<MeshLib::Element*>& elements,
421 std::vector<int>& materials,
422 const std::vector<MeshLib::Node*>& nodes) const
423{
424 std::string line;
425 std::getline(ins, line);
426 std::size_t n_tets;
427 std::size_t n_nodes_per_tet;
428 bool region_attributes;
429
430 while (!ins.fail())
431 {
432 BaseLib::simplify(line);
433 if (line.empty() || line.compare(0, 1, "#") == 0)
434 {
435 // this line is a comment - skip
436 std::getline(ins, line);
437 continue;
438 }
439
440 // read header line
441 bool header_okay = parseElementsFileHeader(
442 line, n_tets, n_nodes_per_tet, region_attributes);
443 if (!header_okay)
444 {
445 return false;
446 }
447 if (!parseElements(ins,
448 elements,
449 materials,
450 nodes,
451 n_tets,
452 n_nodes_per_tet,
453 region_attributes))
454 {
455 return false;
456 }
457 return true;
458 }
459 return false;
460}
461
463 std::size_t& n_tets,
464 std::size_t& n_nodes_per_tet,
465 bool& region_attribute)
466{
467 std::size_t pos_beg;
468 std::size_t pos_end;
469
470 // number of tetrahedras
471 pos_beg = line.find_first_not_of(' ');
472 pos_end = line.find_first_of(' ', pos_beg);
473 if (pos_beg != std::string::npos && pos_end != std::string::npos)
474 {
476 line.substr(pos_beg, pos_end - pos_beg));
477 }
478 else
479 {
480 ERR("TetGenInterface::parseElementsFileHeader(): Could not read number "
481 "of tetrahedra specified in header.");
482 return false;
483 }
484 // nodes per tet - either 4 or 10
485 pos_beg = line.find_first_not_of(" \t", pos_end);
486 pos_end = line.find_first_of(" \t", pos_beg);
487 n_nodes_per_tet = BaseLib::str2number<std::size_t>(
488 line.substr(pos_beg, pos_end - pos_beg));
489 // region attribute at tetrahedra?
490 pos_beg = line.find_first_not_of(" \t", pos_end);
491 pos_end = line.find_first_of(" \t\n", pos_beg);
492 if (pos_end == std::string::npos)
493 {
494 pos_end = line.size();
495 }
496 region_attribute = line.substr(pos_beg, pos_end - pos_beg) == "1";
497
498 return true;
499}
500
501bool TetGenInterface::parseElements(std::ifstream& ins,
502 std::vector<MeshLib::Element*>& elements,
503 std::vector<int>& materials,
504 const std::vector<MeshLib::Node*>& nodes,
505 std::size_t n_tets,
506 std::size_t n_nodes_per_tet,
507 bool region_attribute) const
508{
509 std::string line;
510 std::vector<std::size_t> ids(n_nodes_per_tet);
511
512 elements.reserve(n_tets);
513 materials.reserve(n_tets);
514
515 const unsigned offset = (_zero_based_idx) ? 0 : 1;
516 for (std::size_t k(0); k < n_tets && !ins.fail(); k++)
517 {
518 std::getline(ins, line);
519 if (ins.fail())
520 {
521 ERR("TetGenInterface::parseElements(): Error reading tetrahedron "
522 "{:d}.",
523 k);
524 return false;
525 }
526
527 std::size_t pos_end = 0;
528 std::size_t pos_beg = line.find_first_not_of(' ', pos_end);
529 pos_end = line.find_first_of(" \n", pos_beg);
530
531 if (line.empty() || pos_beg == pos_end ||
532 line.compare(pos_beg, 1, "#") == 0)
533 {
534 k--;
535 continue;
536 }
537
538 if (pos_beg == std::string::npos || pos_end == std::string::npos)
539 {
540 ERR("TetGenInterface::parseElements(): Error reading id of "
541 "tetrahedron {:d}.",
542 k);
543 return false;
544 }
545
546 // read node ids
547 for (std::size_t i(0); i < n_nodes_per_tet; i++)
548 {
549 pos_beg = line.find_first_not_of(' ', pos_end);
550 pos_end = line.find_first_of(' ', pos_beg);
551 if (pos_end == std::string::npos)
552 {
553 pos_end = line.size();
554 }
555 if (pos_beg != std::string::npos && pos_end != std::string::npos)
556 {
558 line.substr(pos_beg, pos_end - pos_beg)) -
559 offset;
560 }
561 else
562 {
563 ERR("TetGenInterface::parseElements(): Error reading node {:d} "
564 "of tetrahedron {:d}.",
565 i,
566 k);
567 return false;
568 }
569 }
570
571 // read region attribute - this is something like material group
572 int region(0);
573 if (region_attribute)
574 {
575 pos_beg = line.find_first_not_of(' ', pos_end);
576 pos_end = line.find_first_of(' ', pos_beg);
577 if (pos_end == std::string::npos)
578 {
579 pos_end = line.size();
580 }
581 if (pos_beg != std::string::npos && pos_end != std::string::npos)
582 {
584 line.substr(pos_beg, pos_end - pos_beg));
585 }
586 else
587 {
588 ERR("TetGenInterface::parseElements(): Error reading region "
589 "attribute of tetrahedron {:d}.",
590 k);
591 return false;
592 }
593 }
594 // insert new element into vector
595 auto** tet_nodes = new MeshLib::Node*[4];
596 for (int n = 0; n < 4; n++)
597 {
598 tet_nodes[n] = nodes[ids[n]];
599 }
600 elements.push_back(new MeshLib::Tet(tet_nodes));
601 materials.push_back(region);
602 }
603
604 return true;
605}
606
608 const std::string& file_name,
609 const GeoLib::GEOObjects& geo_objects,
610 const std::string& geo_name,
611 const std::vector<GeoLib::Point>& attribute_points)
612{
613 std::vector<GeoLib::Point*> const* const points =
614 geo_objects.getPointVec(geo_name);
615 std::vector<GeoLib::Surface*> const* const surfaces =
616 geo_objects.getSurfaceVec(geo_name);
617
618 if (points == nullptr)
619 {
620 ERR("Geometry {:s} not found.", geo_name);
621 return false;
622 }
623 if (surfaces == nullptr)
624 {
625 WARN("No surfaces found for geometry {:s}. Writing points only.",
626 geo_name);
627 }
628
629 std::ofstream out(file_name.c_str(), std::ios::out);
630 out.precision(std::numeric_limits<double>::max_digits10);
631 // the points header
632 const std::size_t nPoints(points->size());
633 out << nPoints << " 3\n";
634 // the point list
635 for (std::size_t i = 0; i < nPoints; ++i)
636 {
637 out << i << " " << (*(*points)[i])[0] << " " << (*(*points)[i])[1]
638 << " " << (*(*points)[i])[2] << "\n";
639 }
640 // the surfaces header
641 const std::size_t nSurfaces = (surfaces) ? surfaces->size() : 0;
642 std::size_t nTotalTriangles(0);
643 for (std::size_t i = 0; i < nSurfaces; ++i)
644 {
645 nTotalTriangles += (*surfaces)[i]->getNumberOfTriangles();
646 }
647 out << nTotalTriangles << " 1\n";
648
649 for (std::size_t i = 0; i < nSurfaces; ++i)
650 {
651 const std::size_t nTriangles((*surfaces)[i]->getNumberOfTriangles());
652 const std::size_t marker(i + 1); // must NOT be 0!
653 // the poly list
654 for (std::size_t j = 0; j < nTriangles; ++j)
655 {
656 const GeoLib::Triangle& tri = *(*(*surfaces)[i])[j];
657 out << "3 " << tri[0] << " " << tri[1] << " " << tri[2] << " "
658 << marker << "\n";
659 }
660 }
661 out << "0\n"; // the polygon holes list
662 // the region attributes list
663 if (attribute_points.empty())
664 {
665 out << "0\n";
666 }
667 else
668 {
669 const std::size_t nAttributePoints(attribute_points.size());
670 out << nAttributePoints << "\n";
671 for (std::size_t i = 0; i < nAttributePoints; ++i)
672 {
673 out << i + 1 << " " << attribute_points[i][0] << " "
674 << attribute_points[i][1] << " " << attribute_points[i][2]
675 << " " << 10 * attribute_points[i].getID() << "\n";
676 }
677 }
678 INFO(
679 "TetGenInterface::writeTetGenSmesh() - {:d} points and {:d} surfaces "
680 "successfully written.",
681 nPoints,
682 nSurfaces);
683 out.close();
684 return true;
685}
686
688 const std::string& file_name,
689 const MeshLib::Mesh& mesh,
690 std::vector<MeshLib::Node>& attribute_points) const
691{
692 if (mesh.getDimension() == 1)
693 {
694 return false;
695 }
696
697 const std::vector<MeshLib::Node*>& nodes = mesh.getNodes();
698
699 std::ofstream out(file_name.c_str(), std::ios::out);
700 out.precision(std::numeric_limits<double>::max_digits10);
701 // the points header
702 const std::size_t nPoints(nodes.size());
703 out << nPoints << " 3\n";
704 // the point list
705 for (std::size_t i = 0; i < nPoints; ++i)
706 {
707 out << i << " " << (*nodes[i])[0] << " " << (*nodes[i])[1] << " "
708 << (*nodes[i])[2] << "\n";
709 }
710
711 if (mesh.getDimension() == 2)
712 {
713 write2dElements(out, mesh);
714 }
715 else
716 {
717 write3dElements(out, mesh, attribute_points);
718 }
719
720 out << "0\n"; // the polygon holes list
721
722 // the region attributes list
723 if (attribute_points.empty())
724 {
725 out << "0\n";
726 }
727 else
728 {
729 const std::size_t nAttributePoints(attribute_points.size());
730 out << nAttributePoints << "\n";
731 for (std::size_t i = 0; i < nAttributePoints; ++i)
732 {
733 out << i + 1 << " " << attribute_points[i][0] << " "
734 << attribute_points[i][1] << " " << attribute_points[i][2]
735 << " " << 10 * attribute_points[i].getID() << "\n";
736 }
737 }
738
739 INFO(
740 "TetGenInterface::writeTetGenPoly() - {:d} points and {:d} surfaces "
741 "successfully written.",
742 nPoints,
743 mesh.getNumberOfElements());
744 out.close();
745 return true;
746}
747
748void TetGenInterface::write2dElements(std::ofstream& out,
749 const MeshLib::Mesh& mesh) const
750{
751 // the surfaces header
752 auto const& types =
754 std::size_t const n_tri =
755 (types.find(MeshLib::MeshElemType::TRIANGLE) != types.end())
757 : 0;
758 std::size_t const n_quad =
759 (types.find(MeshLib::MeshElemType::QUAD) != types.end())
761 : 0;
762 unsigned const nTotalTriangles = n_tri + 2 * n_quad;
763 out << nTotalTriangles << " 1\n";
764
765 const std::vector<MeshLib::Element*>& elements = mesh.getElements();
766 MeshLib::PropertyVector<int> const* const mat_ids =
767 mesh.getProperties().existsPropertyVector<int>("MaterialIDs")
768 ? mesh.getProperties().getPropertyVector<int>("MaterialIDs")
769 : nullptr;
770 const std::size_t nElements(elements.size());
771 unsigned element_count(0);
772 for (std::size_t i = 0; i < nElements; ++i)
773 {
774 std::string const matId = mat_ids ? std::to_string((*mat_ids)[i]) : "";
775 this->writeElementToFacets(out, *elements[i], element_count, matId);
776 }
777}
778
780 std::ofstream& out,
781 const MeshLib::Mesh& mesh,
782 std::vector<MeshLib::Node>& attribute_points) const
783{
784 const std::vector<MeshLib::Element*>& elements = mesh.getElements();
785 const std::size_t nElements(elements.size());
786 if (!attribute_points.empty())
787 {
788 attribute_points.clear();
789 }
790
791 // get position where number of facets need to be written and figure out
792 // worst case of chars that are needed
793 const std::streamoff before_elems_pos(out.tellp());
794 const unsigned n_spaces(
795 static_cast<unsigned>(std::floor(log(nElements * 8))) + 1);
796 out << std::string(n_spaces, ' ') << " 1\n";
797 auto const* const materialIds =
798 mesh.getProperties().getPropertyVector<int>("MaterialIDs");
799 unsigned element_count(0);
800 for (std::size_t i = 0; i < nElements; ++i)
801 {
802 if (elements[i]->getDimension() < 3)
803 {
804 continue;
805 }
806
807 const unsigned nFaces(elements[i]->getNumberOfNeighbors());
808 std::string const mat_id_str =
809 materialIds ? std::to_string((*materialIds)[i]) : "";
810 for (std::size_t j = 0; j < nFaces; ++j)
811 {
812 MeshLib::Element const* const neighbor(elements[i]->getNeighbor(j));
813
814 if (neighbor && materialIds &&
815 (*materialIds)[i] <= (*materialIds)[neighbor->getID()])
816 {
817 continue;
818 }
819
820 std::unique_ptr<MeshLib::Element const> const face(
821 elements[i]->getFace(j));
822 this->writeElementToFacets(out, *face, element_count, mat_id_str);
823 }
824 if (materialIds)
825 {
826 attribute_points.emplace_back(
827 MeshLib::getCenterOfGravity(*elements[i]).data(),
828 (*materialIds)[i]);
829 }
830 }
831 // add number of facets at correct position and jump back
832 const std::streamoff after_elems_pos(out.tellp());
833 out.seekp(before_elems_pos);
834 out << element_count;
835 out.seekp(after_elems_pos);
836}
837
839 const MeshLib::Element& element,
840 unsigned& element_count,
841 std::string const& matId)
842{
843 element_count++;
845 {
846 out << "3 " << getNodeIndex(element, 0) << " "
847 << getNodeIndex(element, 1) << " " << getNodeIndex(element, 2)
848 << " " << matId << " # " << element_count << "\n";
849 }
850 else if (element.getGeomType() == MeshLib::MeshElemType::QUAD)
851 {
852 out << "3 " << getNodeIndex(element, 0) << " "
853 << getNodeIndex(element, 1) << " " << getNodeIndex(element, 2)
854 << " " << matId << " # " << element_count << "\n";
855 element_count++;
856 out << "3 " << getNodeIndex(element, 0) << " "
857 << getNodeIndex(element, 2) << " " << getNodeIndex(element, 3)
858 << " " << matId << " # " << element_count << "\n";
859 }
860}
861
862} // end namespace FileIO
Definition of the Element class.
Filename manipulation routines.
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:45
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:40
Definition of the MeshInformation class.
Definition of the Mesh class.
Definition of the Node class.
Definition of string helper functions.
Definition of the TetGenInterface class.
Definition of the Tet class.
bool readElementsFromStream(std::ifstream &ins, std::vector< MeshLib::Element * > &elements, std::vector< int > &materials, const std::vector< MeshLib::Node * > &nodes) const
bool _boundary_markers
true if boundary markers are set, false otherwise
MeshLib::Mesh * readTetGenMesh(std::string const &nodes_fname, std::string const &ele_fname)
bool parseNodes(std::ifstream &ins, std::vector< MeshLib::Node * > &nodes, std::size_t n_nodes, std::size_t dim)
static bool parseElementsFileHeader(std::string &line, std::size_t &n_tets, std::size_t &n_nodes_per_tet, bool &region_attribute)
void write3dElements(std::ofstream &out, const MeshLib::Mesh &mesh, std::vector< MeshLib::Node > &attribute_points) const
bool readNodesFromStream(std::ifstream &ins, std::vector< MeshLib::Node * > &nodes)
static bool parseNodesFileHeader(std::string const &line, std::size_t &n_nodes, std::size_t &dim, std::size_t &n_attributes, bool &boundary_markers)
bool _zero_based_idx
the value is true if the indexing is zero based, else false
bool readTetGenGeometry(std::string const &geo_fname, GeoLib::GEOObjects &geo_objects)
void write2dElements(std::ofstream &out, const MeshLib::Mesh &mesh) const
bool parseElements(std::ifstream &ins, std::vector< MeshLib::Element * > &elements, std::vector< int > &materials, const std::vector< MeshLib::Node * > &nodes, std::size_t n_tets, std::size_t n_nodes_per_tet, bool region_attribute) const
static bool writeTetGenSmesh(const std::string &file_name, const GeoLib::GEOObjects &geo_objects, const std::string &geo_name, const std::vector< GeoLib::Point > &attribute_points)
static void writeElementToFacets(std::ofstream &out, const MeshLib::Element &element, unsigned &element_count, std::string const &matId)
bool parseSmeshFacets(std::ifstream &input, std::vector< GeoLib::Surface * > &surfaces, const std::vector< GeoLib::Point * > &points, const std::vector< std::size_t > &pnt_id_map)
std::size_t getNFacets(std::ifstream &input)
Returns the declared number of facets in the poly file.
Container class for geometric objects.
Definition GEOObjects.h:57
const std::vector< Point * > * getPointVec(const std::string &name) const
void addPointVec(std::vector< Point * > &&points, std::string &name, PointVec::NameIdMap &&pnt_id_name_map, double const eps=std::sqrt(std::numeric_limits< double >::epsilon()))
const PointVec * getPointVecObj(const std::string &name) const
const std::vector< Surface * > * getSurfaceVec(const std::string &name) const
Returns the surface vector with the given name as a const.
void addSurfaceVec(std::vector< Surface * > &&sfc, const std::string &name, SurfaceVec::NameIdMap &&sfc_names)
const std::vector< std::size_t > & getIDMap() const
Definition PointVec.h:97
A Surface is represented by Triangles. It consists of a reference to a vector of (pointers to) points...
Definition Surface.h:33
std::map< std::string, std::size_t > NameIdMap
Definition TemplateVec.h:41
Class Triangle consists of a reference to a point vector and a vector that stores the indices in the ...
Definition Triangle.h:27
virtual MeshElemType getGeomType() const =0
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:91
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:108
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:111
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition Mesh.h:90
Properties & getProperties()
Definition Mesh.h:136
std::size_t getNumberOfElements() const
Get the number of elements.
Definition Mesh.h:99
Property manager on mesh items. Class Properties manages scalar, vector or matrix properties....
Definition Properties.h:33
bool existsPropertyVector(std::string_view name) const
Definition Properties.h:93
PropertyVector< T > * createNewPropertyVector(std::string_view name, MeshItemType mesh_item_type, std::size_t n_components=1)
Definition Properties.h:17
PropertyVector< T > const * getPropertyVector(std::string_view name) const
static std::map< MeshLib::MeshElemType, unsigned > getNumberOfElementTypes(const MeshLib::Mesh &mesh)
void simplify(std::string &str)
std::string getFileExtension(const std::string &path)
void cleanupVectorElements(std::vector< T * > &items)
Definition Algorithm.h:256
std::string extractBaseNameWithoutExtension(std::string const &pathname)
std::vector< std::string > splitString(std::string const &str)
T str2number(const std::string &str)
Definition StringTools.h:63
auto constructPointsFromNodes(std::vector< MeshLib::Node * > nodes)
MathLib::Point3d getCenterOfGravity(Element const &element)
Calculates the center of gravity for the mesh element.
Definition Element.cpp:142