From 482a2aec30c8a04ecaed936e74bc312c0c7ed03f Mon Sep 17 00:00:00 2001 From: Bernd Hahnebach Date: Mon, 19 Mar 2018 13:01:37 +0100 Subject: [PATCH] FEM: VTK export, fix vtk mesh builder and vtk result obj --- src/Mod/Fem/App/FemVTKTools.cpp | 135 ++++++++++++++++++++------------ 1 file changed, 84 insertions(+), 51 deletions(-) diff --git a/src/Mod/Fem/App/FemVTKTools.cpp b/src/Mod/Fem/App/FemVTKTools.cpp index 0485fda46a..1fcc1963eb 100644 --- a/src/Mod/Fem/App/FemVTKTools.cpp +++ b/src/Mod/Fem/App/FemVTKTools.cpp @@ -75,16 +75,18 @@ #include #include -#include -#include -#include -#include -#include -#include #include #include #include #include +#include +#include +#include +#include +#include +#include +#include +#include #include "FemVTKTools.h" #include "FemMeshProperty.h" @@ -212,6 +214,7 @@ FemMesh* FemVTKTools::readVTKMesh(const char* filename, FemMesh* mesh) void exportFemMeshFaces(vtkSmartPointer grid, const SMDS_FaceIteratorPtr& aFaceIter) { + Base::Console().Log(" Start: VTK mesh builder faces.\n"); vtkSmartPointer triangleArray = vtkSmartPointer::New(); vtkSmartPointer quadTriangleArray = vtkSmartPointer::New(); @@ -274,57 +277,60 @@ void exportFemMeshFaces(vtkSmartPointer grid, const SMDS_Fa { throw std::runtime_error("Face not yet supported by FreeCAD's VTK mesh builder\n"); } - } - if(triangleArray->GetNumberOfCells()>0) - grid->SetCells(VTK_TRIANGLE, triangleArray); + } + if(triangleArray->GetNumberOfCells()>0) + grid->SetCells(VTK_TRIANGLE, triangleArray); - if(quadArray->GetNumberOfCells()>0) - grid->SetCells(VTK_QUAD, quadArray); + if(quadArray->GetNumberOfCells()>0) + grid->SetCells(VTK_QUAD, quadArray); - if(quadTriangleArray->GetNumberOfCells()>0) - grid->SetCells(VTK_QUADRATIC_TRIANGLE, quadTriangleArray); + if(quadTriangleArray->GetNumberOfCells()>0) + grid->SetCells(VTK_QUADRATIC_TRIANGLE, quadTriangleArray); - if(quadQuadArray->GetNumberOfCells()>0) - grid->SetCells(VTK_QUADRATIC_QUAD, quadQuadArray); + if(quadQuadArray->GetNumberOfCells()>0) + grid->SetCells(VTK_QUADRATIC_QUAD, quadQuadArray); + Base::Console().Log(" End: VTK mesh builder faces.\n"); } void exportFemMeshCells(vtkSmartPointer grid, const SMDS_VolumeIteratorPtr& aVolIter) { - // add common CFD cells like hex, wedge, prism, in addition to tetra + Base::Console().Log(" Start: VTK mesh builder volumes.\n"); + vtkSmartPointer tetraArray = vtkSmartPointer::New(); - vtkSmartPointer hexaArray = vtkSmartPointer::New(); - vtkSmartPointer wedgeArray = vtkSmartPointer::New(); vtkSmartPointer pyramidArray = vtkSmartPointer::New(); - // quadratic elemnts with 13 and 15 nodes are not added yet + vtkSmartPointer wedgeArray = vtkSmartPointer::New(); + vtkSmartPointer hexaArray = vtkSmartPointer::New(); vtkSmartPointer quadTetraArray = vtkSmartPointer::New(); + vtkSmartPointer quadPyramidArray = vtkSmartPointer::New(); + vtkSmartPointer quadWedgeArray = vtkSmartPointer::New(); vtkSmartPointer quadHexaArray = vtkSmartPointer::New(); for (;aVolIter->more();) { const SMDS_MeshVolume* aVol = aVolIter->next(); - if (aVol->NbNodes() == 4) { // tetrahedra - vtkSmartPointer tetra = vtkSmartPointer::New(); - tetra->GetPointIds()->SetId(0, aVol->GetNode(0)->GetID()-1); - tetra->GetPointIds()->SetId(1, aVol->GetNode(1)->GetID()-1); - tetra->GetPointIds()->SetId(2, aVol->GetNode(2)->GetID()-1); - tetra->GetPointIds()->SetId(3, aVol->GetNode(3)->GetID()-1); - - tetraArray->InsertNextCell(tetra); + if (aVol->NbNodes() == 4) { // tetra4 + Base::Console().Log(" Volume tetra4\n"); + vtkSmartPointer cell = vtkSmartPointer::New(); + cell->GetPointIds()->SetId(0, aVol->GetNode(0)->GetID()-1); + cell->GetPointIds()->SetId(1, aVol->GetNode(1)->GetID()-1); + cell->GetPointIds()->SetId(2, aVol->GetNode(2)->GetID()-1); + cell->GetPointIds()->SetId(3, aVol->GetNode(3)->GetID()-1); + tetraArray->InsertNextCell(cell); } - else if (aVol->NbNodes() == 5) { // common cell types for CFD - vtkSmartPointer cell= vtkSmartPointer::New(); + else if (aVol->NbNodes() == 5) { // pyra5 + Base::Console().Log(" Volume pyra5\n"); + vtkSmartPointer cell = vtkSmartPointer::New(); cell->GetPointIds()->SetId(0, aVol->GetNode(0)->GetID()-1); cell->GetPointIds()->SetId(1, aVol->GetNode(1)->GetID()-1); cell->GetPointIds()->SetId(2, aVol->GetNode(2)->GetID()-1); cell->GetPointIds()->SetId(3, aVol->GetNode(3)->GetID()-1); cell->GetPointIds()->SetId(4, aVol->GetNode(4)->GetID()-1); - pyramidArray->InsertNextCell(cell); } - /* // penta6 and hexa8 as well has hexa20 make FreeCAD crash - else if (aVol->NbNodes() == 6) { + else if (aVol->NbNodes() == 6) { // penta6 + Base::Console().Log(" Volume penta6\n"); vtkSmartPointer cell = vtkSmartPointer::New(); cell->GetPointIds()->SetId(0, aVol->GetNode(0)->GetID()-1); cell->GetPointIds()->SetId(1, aVol->GetNode(1)->GetID()-1); @@ -332,11 +338,11 @@ void exportFemMeshCells(vtkSmartPointer grid, const SMDS_Vo cell->GetPointIds()->SetId(3, aVol->GetNode(3)->GetID()-1); cell->GetPointIds()->SetId(4, aVol->GetNode(4)->GetID()-1); cell->GetPointIds()->SetId(5, aVol->GetNode(5)->GetID()-1); - wedgeArray->InsertNextCell(cell); } - else if (aVol->NbNodes() == 8) { - vtkSmartPointer cell= vtkSmartPointer::New(); + else if (aVol->NbNodes() == 8) { // hexa8 + Base::Console().Log(" Volume hexa8\n"); + vtkSmartPointer cell = vtkSmartPointer::New(); cell->GetPointIds()->SetId(0, aVol->GetNode(0)->GetID()-1); cell->GetPointIds()->SetId(1, aVol->GetNode(1)->GetID()-1); cell->GetPointIds()->SetId(2, aVol->GetNode(2)->GetID()-1); @@ -345,28 +351,42 @@ void exportFemMeshCells(vtkSmartPointer grid, const SMDS_Vo cell->GetPointIds()->SetId(5, aVol->GetNode(5)->GetID()-1); cell->GetPointIds()->SetId(6, aVol->GetNode(6)->GetID()-1); cell->GetPointIds()->SetId(7, aVol->GetNode(7)->GetID()-1); - hexaArray->InsertNextCell(cell); } - */ - else if (aVol->NbNodes() == 10) { // quadratic tetrahedra - + else if (aVol->NbNodes() == 10) { // tetra10 + Base::Console().Log(" Volume tetra10\n"); vtkSmartPointer tetra = vtkSmartPointer::New(); for(int i=0; i<10; i++){ tetra->GetPointIds()->SetId(i, aVol->GetNode(i)->GetID()-1); } - quadTetraArray->InsertNextCell(tetra); } - /* - else if (aVol->NbNodes() == 20) { // not tested, no sure about vertex sequence - vtkSmartPointer cell= vtkSmartPointer::New(); + + else if (aVol->NbNodes() == 13) { // pyra13 + Base::Console().Log(" Volume pyra13\n"); + vtkSmartPointer cell = vtkSmartPointer::New(); + for(int i=0; i<13; i++){ + cell->GetPointIds()->SetId(i, aVol->GetNode(i)->GetID()-1); + // Base::Console().Log("node ids: %i\n", aVol->GetNode(i)->GetID()-1); + } + quadPyramidArray->InsertNextCell(cell); + } + else if (aVol->NbNodes() == 15) { // penta15 + Base::Console().Log(" Volume penta15\n"); + vtkSmartPointer cell = vtkSmartPointer::New(); + for(int i=0; i<15; i++){ + cell->GetPointIds()->SetId(i, aVol->GetNode(i)->GetID()-1); + } + quadWedgeArray->InsertNextCell(cell); + } + else if (aVol->NbNodes() == 20) { // hexa20 + Base::Console().Log(" Volume hexa20\n"); + vtkSmartPointer cell = vtkSmartPointer::New(); for(int i=0; i<20; i++){ cell->GetPointIds()->SetId(i, aVol->GetNode(i)->GetID()-1); } quadHexaArray->InsertNextCell(cell); } - */ else { throw std::runtime_error("Volume not yet supported by FreeCAD's VTK mesh builder\n"); } @@ -387,14 +407,22 @@ void exportFemMeshCells(vtkSmartPointer grid, const SMDS_Vo if(quadTetraArray->GetNumberOfCells()>0) grid->SetCells(VTK_QUADRATIC_TETRA, quadTetraArray); + if(quadPyramidArray->GetNumberOfCells()>0) + grid->SetCells(VTK_QUADRATIC_PYRAMID, quadPyramidArray); + + if(quadWedgeArray->GetNumberOfCells()>0) + grid->SetCells(VTK_QUADRATIC_WEDGE, quadWedgeArray); + if(quadHexaArray->GetNumberOfCells()>0) grid->SetCells(VTK_QUADRATIC_HEXAHEDRON, quadHexaArray); + Base::Console().Log(" End: VTK mesh builder volumes.\n"); } void FemVTKTools::exportVTKMesh(const FemMesh* mesh, vtkSmartPointer grid, float scale) { + Base::Console().Message("Start: VTK mesh builder ======================\n"); SMESH_Mesh* smesh = const_cast(mesh->getSMesh()); SMESHDS_Mesh* meshDS = smesh->GetMeshDS(); const SMDS_MeshInfo& info = meshDS->GetMeshInfo(); @@ -410,20 +438,21 @@ void FemVTKTools::exportVTKMesh(const FemMesh* mesh, vtkSmartPointerSetPoint(node->GetID()-1, coords); } grid->SetPoints(points); - //start with 2d elements + // faces SMDS_FaceIteratorPtr aFaceIter = meshDS->facesIterator(); exportFemMeshFaces(grid, aFaceIter); - //3D volume elements + // volumes SMDS_VolumeIteratorPtr aVolIter = meshDS->volumesIterator(); exportFemMeshCells(grid, aVolIter); + Base::Console().Message("End: VTK mesh builder ======================\n"); } void FemVTKTools::writeVTKMesh(const char* filename, const FemMesh* mesh) { Base::TimeInfo Start; - Base::Console().Log("Start: write FemMesh from VTK unstructuredGrid ======================\n"); + Base::Console().Message("Start: write FemMesh from VTK unstructuredGrid ======================\n"); Base::FileInfo f(filename); vtkSmartPointer grid = vtkSmartPointer::New(); @@ -439,7 +468,7 @@ void FemVTKTools::writeVTKMesh(const char* filename, const FemMesh* mesh) Base::Console().Error("file name extension is not supported to write VTK\n"); } - Base::Console().Log(" %f: Done \n",Base::TimeInfo::diffTimeF(Start, Base::TimeInfo())); + Base::Console().Message(" %f: Done \n",Base::TimeInfo::diffTimeF(Start, Base::TimeInfo())); } @@ -589,13 +618,15 @@ void FemVTKTools::writeResult(const char* filename, const App::DocumentObject* r return; } + float scale = 1.0; + /* auto hGrp = App::GetApplication().GetParameterGroupByPath("User parameter:BaseApp/Preferences/Units"); int unitSchema = hGrp->GetInt("UserSchema",0); - float scale = 1.0; if(unitSchema == 0) // standard mm { - scale = 0.001; // convert from mm in FreeCAD to SI length in result file + //scale = 0.001; // convert from mm in FreeCAD to SI length in result file } + */ Base::TimeInfo Start; Base::Console().Message("Start: write FemResult or CfdResult to VTK unstructuredGrid dataset =======\n"); @@ -606,6 +637,8 @@ void FemVTKTools::writeResult(const char* filename, const App::DocumentObject* r const FemMesh& fmesh = static_cast(mesh->getPropertyByName("FemMesh"))->getValue(); FemVTKTools::exportVTKMesh(&fmesh, grid, scale); + Base::Console().Message(" %f: vtk mesh builder finisched\n",Base::TimeInfo::diffTimeF(Start, Base::TimeInfo())); + if(res->getPropertyByName("Velocity")){ // consider better way to detect result type, res->Type == "CfdResult" FemVTKTools::exportFluidicResult(res, grid); } @@ -628,7 +661,7 @@ void FemVTKTools::writeResult(const char* filename, const App::DocumentObject* r Base::Console().Error("file name extension is not supported to write VTK\n"); } - Base::Console().Message(" %f: result writing is done \n",Base::TimeInfo::diffTimeF(Start, Base::TimeInfo())); + Base::Console().Message(" %f: writing result object to vtk finished\n",Base::TimeInfo::diffTimeF(Start, Base::TimeInfo())); } // it is an internal utility func to avoid code duplication