diff --git a/FileIO/TetGenInterface.cpp b/FileIO/TetGenInterface.cpp index 7eb739682956001949c5c481c715d77e5b9b4bf6..0fe6af93616f3cf59223c945aaa96b0deb9ee016 100644 --- a/FileIO/TetGenInterface.cpp +++ b/FileIO/TetGenInterface.cpp @@ -58,9 +58,9 @@ bool TetGenInterface::readTetGenGeometry (std::string const& geo_fname, return false; } std::string ext (BaseLib::getFileExtension(geo_fname)); - if (ext.compare("poly") && ext.compare("smesh")) + if (ext.compare("smesh")) { - ERR ("TetGenInterface::readTetGenPoly() - unknown file type (only *.poly or *.smesh supported)."); + ERR ("TetGenInterface::readTetGenPoly() - unknown file type (only *.smesh are supported)."); return false; } @@ -85,12 +85,7 @@ bool TetGenInterface::readTetGenGeometry (std::string const& geo_fname, const std::vector<std::size_t> &id_map (geo_objects.getPointVecObj(geo_name)->getIDMap()); std::vector<GeoLib::Surface*> *surfaces = new std::vector<GeoLib::Surface*>; - bool failed (true); - if (ext.compare("poly") == 0) - failed = !parsePolyFacets(poly_stream, *surfaces, *points, id_map); - else - failed = !parseSmeshFacets(poly_stream, *surfaces, *points, id_map); - if (failed) + if (!parseSmeshFacets(poly_stream, *surfaces, *points, id_map)) { // remove surfaces read until now but keep the points for (std::size_t k=0; k<surfaces->size(); k++) @@ -130,85 +125,6 @@ std::size_t TetGenInterface::getNFacets(std::ifstream &input) return false; } -bool TetGenInterface::parsePolyFacets(std::ifstream &input, - std::vector<GeoLib::Surface*> &surfaces, - const std::vector<GeoLib::Point*> &points, - const std::vector<std::size_t> &pnt_id_map) -{ - const std::size_t nFacets (this->getNFacets(input)); - std::size_t nMultPolys (0); - std::string line; - surfaces.reserve(nFacets); - std::list<std::string>::const_iterator it; - - const unsigned offset = (_zero_based_idx) ? 0 : 1; - for (std::size_t k(0); k<nFacets && !input.fail(); k++) - { - getline (input, line); - if (input.fail()) - { - ERR("TetGenInterface::parseFacets(): Error reading facet %d.", k); - return false; - } - - BaseLib::simplify(line); - if (line.empty() || line.compare(0,1,"#") == 0) - { - k--; - continue; - } - - // read facets - const std::list<std::string> poly_def_fields = BaseLib::splitString(line, ' '); - it = poly_def_fields.begin(); - const std::size_t nPolys = BaseLib::str2number<std::size_t>(*it); - const std::size_t nPolyHoles = (poly_def_fields.size()>1) ? BaseLib::str2number<std::size_t>(*(++it)) : 0; - // here this line also potentially includes a boundary marker which we ignore for now - nMultPolys += (nPolys-1); - - // read polys - for (std::size_t i(0); i<nPolys && !input.fail(); ++i) - { - getline (input, line); - BaseLib::simplify(line); - if (line.empty() || line.compare(0,1,"#") == 0) - { - i--; - continue; - } - - const std::list<std::string> point_fields = BaseLib::splitString(line, ' '); - it = point_fields.begin(); - const std::size_t nPoints = BaseLib::str2number<std::size_t>(*it); - if (point_fields.size() > nPoints) - { - GeoLib::Polyline polyline(points); - for (std::size_t j(0); j<nPoints; ++j) - polyline.addPoint(pnt_id_map[BaseLib::str2number<std::size_t>(*(++it))-offset]); - - polyline.closePolyline(); - surfaces.push_back(GeoLib::Surface::createSurface(polyline)); - } - else - { - ERR("TetGenInterface::parseFacets(): Error reading points for polygon %d of facet %d.", i, k); - return false; - } - } - for (std::size_t j(0); j<nPolyHoles && !input.fail(); ++j) - getline(input, line); - // Here are points defined which are located in holes within the surface. We ignore these as they are not part of the actual geometry. - } - // here the poly-file potentially defines a number of points to mark holes within the volumes defined by the facets, these are ignored for now - // here the poly-file potentially defines a number of region attributes, these are ignored for now - - if (surfaces.size() == nFacets+nMultPolys) - return true; - - ERR ("TetGenInterface::parseFacets(): Number of expected surfaces (%d) does not match number of found surfaces (%d).", nFacets+nMultPolys, surfaces.size()); - return false; -} - bool TetGenInterface::parseSmeshFacets(std::ifstream &input, std::vector<GeoLib::Surface*> &surfaces, const std::vector<GeoLib::Point*> &points, @@ -640,13 +556,12 @@ bool TetGenInterface::writeTetGenSmesh(const std::string &file_name, bool TetGenInterface::writeTetGenSmesh(const std::string &file_name, const MeshLib::Mesh &mesh, - const std::vector<MeshLib::Node> &attribute_points) const + std::vector<MeshLib::Node> &attribute_points) const { - if (mesh.getDimension() != 2) + if (mesh.getDimension() == 1) return false; const std::vector<MeshLib::Node*> &nodes = mesh.getNodes(); - const std::vector<MeshLib::Element*> &elements = mesh.getElements(); std::ofstream out( file_name.c_str(), std::ios::out ); out.precision(std::numeric_limits<double>::digits10); @@ -657,25 +572,11 @@ bool TetGenInterface::writeTetGenSmesh(const std::string &file_name, for (std::size_t i=0; i<nPoints; ++i) out << i << " " << (*nodes[i])[0] << " " << (*nodes[i])[1] << " " << (*nodes[i])[2] << "\n"; - // the surfaces header - const std::array<unsigned,7> types = MeshInformation::getNumberOfElementTypes(mesh); - const unsigned nTotalTriangles (types[1] + (2*types[2])); - out << nTotalTriangles << " 1\n"; + if (mesh.getDimension() == 2) + write2dElements(out, mesh, attribute_points); + else + write3dElements(out, mesh, attribute_points); - const std::size_t nElements (elements.size()); - unsigned count(0); - for (std::size_t i=0; i<nElements; ++i) - { - count++; - if (elements[i]->getGeomType() == MeshElemType::TRIANGLE) - out << "3 " << elements[i]->getNodeIndex(0) << " " << elements[i]->getNodeIndex(1) << " " << elements[i]->getNodeIndex(2) << " " << elements[i]->getValue() << " # " << count << "\n"; - else if (elements[i]->getGeomType() == MeshElemType::QUAD) - { - out << "3 " << elements[i]->getNodeIndex(0) << " " << elements[i]->getNodeIndex(1) << " " << elements[i]->getNodeIndex(2) << " " << elements[i]->getValue() << " # " << count << "\n"; - count++; - out << "3 " << elements[i]->getNodeIndex(0) << " " << elements[i]->getNodeIndex(2) << " " << elements[i]->getNodeIndex(3) << " " << elements[i]->getValue() << " # " << count << "\n"; - } - } out << "0\n"; // the polygon holes list // the region attributes list @@ -688,9 +589,89 @@ bool TetGenInterface::writeTetGenSmesh(const std::string &file_name, for (std::size_t i=0; i<nAttributePoints; ++i) out << i+1 << " " << attribute_points[i][0] << " " << attribute_points[i][1] << " " << attribute_points[i][2] << " " << 10*attribute_points[i].getID() << "\n"; } - INFO ("TetGenInterface::writeTetGenPoly() - %d points and %d surfaces successfully written.", nPoints, nElements); + + INFO ("TetGenInterface::writeTetGenPoly() - %d points and %d surfaces successfully written.", nPoints, mesh.getNElements()); out.close(); return true; } +void TetGenInterface::write2dElements(std::ofstream &out, + const MeshLib::Mesh &mesh, + std::vector<MeshLib::Node> &attribute_points) const +{ + // the surfaces header + const std::array<unsigned,7> types = MeshInformation::getNumberOfElementTypes(mesh); + const unsigned nTotalTriangles (types[1] + (2*types[2])); + out << nTotalTriangles << " 1\n"; + + const std::vector<MeshLib::Element*> &elements = mesh.getElements(); + const std::size_t nElements (elements.size()); + unsigned element_count(0); + for (std::size_t i=0; i<nElements; ++i) + this->writeElementToFacets(out, *elements[i], element_count); +} + +void TetGenInterface::write3dElements(std::ofstream &out, + const MeshLib::Mesh &mesh, + std::vector<MeshLib::Node> &attribute_points) const +{ + const std::vector<MeshLib::Element*> &elements = mesh.getElements(); + const std::size_t nElements (elements.size()); + if (!attribute_points.empty()) + attribute_points.clear(); + + // get position where number of facets need to be written and figure out worst case of chars that are needed + const std::streamoff before_elems_pos (out.tellp()); + const unsigned n_spaces (static_cast<unsigned>(floor(log(nElements*8))) + 1); + out << std::string(n_spaces, ' ') << "\n"; + + unsigned element_count(0); + for (std::size_t i=0; i<nElements; ++i) + { + if (elements[i]->getDimension() < 3) + continue; + + const unsigned nFaces (elements[i]->getNNeighbors()); + for (std::size_t j=0; j<nFaces; ++j) + { + MeshLib::Element const*const neighbor ( elements[i]->getNeighbor(j) ); + + if (neighbor) + { + if (elements[i]->getValue() > neighbor->getValue()) + { + MeshLib::Element const*const face (elements[i]->getFace(j)); + this->writeElementToFacets(out, *face, element_count); + delete face; + } + } + else + { + MeshLib::Element const*const face (elements[i]->getFace(j)); + this->writeElementToFacets(out, *face, element_count); + delete face; + } + } + attribute_points.push_back(MeshLib::Node(elements[i]->getCenterOfGravity().getCoords(), elements[i]->getValue())); + } + // add number of facets at correct position and jump back + const std::streamoff after_elems_pos (out.tellp()); + out.seekp(before_elems_pos); + out << element_count; + out.seekp(after_elems_pos); +} + +void TetGenInterface::writeElementToFacets(std::ofstream &out, const MeshLib::Element &element, unsigned &element_count) const +{ + element_count++; + if (element.getGeomType() == MeshElemType::TRIANGLE) + out << "3 " << element.getNodeIndex(0) << " " << element.getNodeIndex(1) << " " << element.getNodeIndex(2) << " " << element.getValue() << " # " << element_count << "\n"; + else if (element.getGeomType() == MeshElemType::QUAD) + { + out << "3 " << element.getNodeIndex(0) << " " << element.getNodeIndex(1) << " " << element.getNodeIndex(2) << " " << element.getValue() << " # " << element_count << "\n"; + element_count++; + out << "3 " << element.getNodeIndex(0) << " " << element.getNodeIndex(2) << " " << element.getNodeIndex(3) << " " << element.getValue() << " # " << element_count << "\n"; + } +} + } // end namespace FileIO diff --git a/FileIO/TetGenInterface.h b/FileIO/TetGenInterface.h index 2e5a38555f8f9600e4178c0d323e49d222e7255f..8ed6a9b5e13a930abfd14423e7d05ec576d497bd 100644 --- a/FileIO/TetGenInterface.h +++ b/FileIO/TetGenInterface.h @@ -82,25 +82,12 @@ public: */ bool writeTetGenSmesh(const std::string &file_name, const MeshLib::Mesh &mesh, - const std::vector<MeshLib::Node> &attribute_points) const; + std::vector<MeshLib::Node> &attribute_points) const; private: /// Returns the declared number of facets in the poly file. std::size_t getNFacets(std::ifstream &input); - /** - * Method parses the lines reading the facets from TetGen poly file - * @param input the input stream (input) - * @param surfaces the vector of surfaces to be filled (output) - * @param points the point vector needed for creating surfaces (input) - * @param pnt_id_map the id map to compensate for possibly changed point ids after adding the point vector to GEOObjects - * @return true, if the facets have been read correctly, false if the method detects an error - */ - bool parsePolyFacets(std::ifstream &input, - std::vector<GeoLib::Surface*> &surfaces, - const std::vector<GeoLib::Point*> &points, - const std::vector<std::size_t> &pnt_id_map); - /** * Method parses the lines reading the facets from TetGen smesh file * @param input the input stream (input) @@ -190,6 +177,30 @@ private: std::size_t n_nodes_per_tet, bool region_attribute) const; + /** + * Writes the elements from a 2D mesh to a TetGen smesh-file. + * @param out the output stream the information is written to. + * @param mesh mesh containing the subsurface boundary representation used for meshing. + * @param attribute_points attribute points containing material IDs (if the vector is empty no attributes are written). + * @return returns true on success and false otherwise. + */ + void write2dElements(std::ofstream &out, + const MeshLib::Mesh &mesh, + std::vector<MeshLib::Node> &attribute_points) const; + + /** + * Writes the elements from a 3D mesh to a TetGen smesh-file. + * @param out the output stream the information is written to. + * @param mesh the 3D mesh. + * @param attribute_points attribute points containing material IDs (emptied when called and then filled with correct values). + * @return returns true on success and false otherwise. + */ + void write3dElements(std::ofstream &out, + const MeshLib::Mesh &mesh, + std::vector<MeshLib::Node> &attribute_points) const; + + /// Writes facet information from a 2D element to the stream and increments the total element count accordingly + void writeElementToFacets(std::ofstream &out, const MeshLib::Element &element, unsigned &element_count) const; /// the value is true if the indexing is zero based, else false bool _zero_based_idx; diff --git a/Gui/DataView/MeshLayerEditDialog.cpp b/Gui/DataView/MeshLayerEditDialog.cpp index fc2011ff510686da324a4573323e5f5cb493d651..24d8f4c17c9b37af76436c1c5eecd2870ab445a9 100644 --- a/Gui/DataView/MeshLayerEditDialog.cpp +++ b/Gui/DataView/MeshLayerEditDialog.cpp @@ -214,29 +214,44 @@ MeshLib::Mesh* MeshLayerEditDialog::createPrismMesh() MeshLib::Mesh* MeshLayerEditDialog::createTetMesh() { + QSettings settings; + QString filename = QFileDialog::getSaveFileName(this, "Write TetGen input file to", + settings.value("lastOpenedTetgenFileDirectory").toString(), + "TetGen Geometry (*.smesh)"); + if (filename.isEmpty()) + return nullptr; + const unsigned nLayers = _layerEdit->text().toInt(); - std::vector<std::string> raster_paths(nLayers+1); - for (unsigned i=0; i<=nLayers; ++i) - raster_paths[i] = this->_edits[i+1]->text().toStdString(); - LayeredVolume lv; - lv.createGeoVolumes(*_msh, raster_paths); + MeshLib::Mesh* tg_mesh (nullptr); + if (_use_rasters) + { + std::vector<std::string> raster_paths(nLayers+1); + for (unsigned i=0; i<=nLayers; ++i) + raster_paths[i] = this->_edits[i+1]->text().toStdString(); + LayeredVolume lv; + lv.createGeoVolumes(*_msh, raster_paths); - MeshLib::Mesh* tg_mesh (lv.getMesh()); + tg_mesh = lv.getMesh(); - QString file_path(""); - if (tg_mesh) - { - QSettings settings; - QString filename = QFileDialog::getSaveFileName(this, "Write TetGen input file to", - settings.value("lastOpenedTetgenFileDirectory").toString(), - "TetGen Geometry (*.smesh)"); - if (!filename.isEmpty()) + QString file_path(""); + if (tg_mesh) { std::vector<MeshLib::Node> tg_attr (lv.getAttributePoints()); FileIO::TetGenInterface tetgen_interface; tetgen_interface.writeTetGenSmesh(filename.toStdString(), *tg_mesh, tg_attr); } } + else + { + std::vector<float> layer_thickness; + for (unsigned i=0; i<nLayers; ++i) + layer_thickness.push_back(this->_edits[i]->text().toFloat()); + tg_mesh = MeshLayerMapper::CreateLayers(*_msh, layer_thickness); + std::vector<MeshLib::Node> tg_attr; + FileIO::TetGenInterface tetgen_interface; + tetgen_interface.writeTetGenSmesh(filename.toStdString(), *tg_mesh, tg_attr); + } + return tg_mesh; } diff --git a/Gui/DataView/MshView.cpp b/Gui/DataView/MshView.cpp index 71862bafaab2f0eaab3f4fe754b805b16ac4e352..cd9db9cfd9dfd410a9217955aa2bd06ed9b8880b 100644 --- a/Gui/DataView/MshView.cpp +++ b/Gui/DataView/MshView.cpp @@ -22,6 +22,7 @@ #include <QSettings> #include "Mesh.h" +#include "Node.h" #include "MeshLayerEditDialog.h" #include "MeshValueEditDialog.h" #include "MshItem.h" @@ -39,6 +40,7 @@ #include "XmlIO/Boost/BoostVtuInterface.h" #include "Writer.h" // necessary to avoid Linker Error in Windows #include "SHPInterface.h" +#include "TetGenInterface.h" MshView::MshView( QWidget* parent /*= 0*/ ) : QTreeView(parent) @@ -120,11 +122,15 @@ void MshView::contextMenuEvent( QContextMenuEvent* event ) QAction* editMeshAction = menu.addAction("Edit mesh..."); QAction* editValuesAction = menu.addAction("Edit material groups..."); QAction* meshQualityAction = menu.addAction("Calculate element quality..."); - QAction* surfaceMeshAction (NULL); + QAction* surfaceMeshAction (nullptr); + QAction* tetgenExportAction (nullptr); if (mesh_dim==3) - surfaceMeshAction = menu.addAction("Extract surface"); - QAction* mesh2geoAction (NULL); - QAction* shapeExportAction (NULL); + { + surfaceMeshAction = menu.addAction("Extract surface"); + tetgenExportAction = menu.addAction("Export to TetGen..."); + } + QAction* mesh2geoAction (nullptr); + QAction* shapeExportAction (nullptr); if (mesh_dim==2) { mesh2geoAction = menu.addAction("Convert to geometry"); @@ -136,17 +142,20 @@ void MshView::contextMenuEvent( QContextMenuEvent* event ) QAction* addDirectAction = direct_cond_menu.addAction("Add..."); QAction* loadDirectAction = direct_cond_menu.addAction("Load..."); //menu.addSeparator(); - connect(editMeshAction, SIGNAL(triggered()), this, SLOT(openMeshEditDialog())); - connect(editValuesAction, SIGNAL(triggered()), this, SLOT(openValuesEditDialog())); - connect(meshQualityAction, SIGNAL(triggered()), this, SLOT(checkMeshQuality())); + connect(editMeshAction, SIGNAL(triggered()), this, SLOT(openMeshEditDialog())); + connect(editValuesAction, SIGNAL(triggered()), this, SLOT(openValuesEditDialog())); + connect(meshQualityAction, SIGNAL(triggered()), this, SLOT(checkMeshQuality())); if (mesh_dim==3) - connect(surfaceMeshAction, SIGNAL(triggered()), this, SLOT(extractSurfaceMesh())); - connect(addDirectAction, SIGNAL(triggered()), this, SLOT(addDIRECTSourceTerms())); - connect(loadDirectAction, SIGNAL(triggered()), this, SLOT(loadDIRECTSourceTerms())); + { + connect(surfaceMeshAction, SIGNAL(triggered()), this, SLOT(extractSurfaceMesh())); + connect(tetgenExportAction, SIGNAL(triggered()), this, SLOT(exportToTetGen())); + } + connect(addDirectAction, SIGNAL(triggered()), this, SLOT(addDIRECTSourceTerms())); + connect(loadDirectAction, SIGNAL(triggered()), this, SLOT(loadDIRECTSourceTerms())); if (mesh_dim==2) { - connect(mesh2geoAction, SIGNAL(triggered()), this, SLOT(convertMeshToGeometry())); - connect(shapeExportAction, SIGNAL(triggered()), this, SLOT(exportToShapefile())); + connect(mesh2geoAction, SIGNAL(triggered()), this, SLOT(convertMeshToGeometry())); + connect(shapeExportAction, SIGNAL(triggered()), this, SLOT(exportToShapefile())); } menu.exec(event->globalPos()); } @@ -212,6 +221,26 @@ void MshView::exportToShapefile() const OGSError::box("Error exporting mesh\n to shapefile"); } +void MshView::exportToTetGen() +{ + QModelIndex index = this->selectionModel()->currentIndex(); + + if (!index.isValid()) + return; + + const MeshLib::Mesh* mesh = static_cast<MshModel*>(this->model())->getMesh(index); + QSettings settings; + QString filename = QFileDialog::getSaveFileName(this, "Write TetGen input file to", + settings.value("lastOpenedTetgenFileDirectory").toString(), + "TetGen Geometry (*.smesh)"); + if (!filename.isEmpty()) + { + FileIO::TetGenInterface tg; + std::vector<MeshLib::Node> attr; + tg.writeTetGenSmesh(filename.toStdString(), *mesh, attr); + } +} + int MshView::writeToFile() const { QModelIndex index = this->selectionModel()->currentIndex(); diff --git a/Gui/DataView/MshView.h b/Gui/DataView/MshView.h index 124037ab295ea0cb1b401bdb59b10470db73fcc4..4ff9c3e700db31ab8a0c7083b9289b443023a90e 100644 --- a/Gui/DataView/MshView.h +++ b/Gui/DataView/MshView.h @@ -73,6 +73,8 @@ private slots: void extractSurfaceMesh(); + void exportToTetGen(); + void loadDIRECTSourceTerms(); void convertMeshToGeometry();