kindred/create
1
1
Fork 0
Code Issues 51 Pull Requests 1 Actions 1 Packages Projects 9 Releases 2 Wiki Activity

FEM: VTK export, fix vtk mesh builder and vtk result obj

Browse Source
This commit is contained in:
Bernd Hahnebach
2018-03-19 13:01:37 +01:00
committed by wmayer
parent 9bf2fe8e38
commit 482a2aec30
1 changed files with 84 additions and 51 deletions
Download Patch File Download Diff File Expand all files Collapse all files

135
src/Mod/Fem/App/FemVTKTools.cpp
View File

@@ -75,16 +75,18 @@
#include <vtkIdList.h>
#include <vtkCellTypes.h>
#include <vtkTetra.h>
#include <vtkHexahedron.h>
#include <vtkWedge.h>
#include <vtkPyramid.h>
#include <vtkQuadraticTetra.h>
#include <vtkQuadraticHexahedron.h>
#include <vtkTriangle.h>
#include <vtkQuad.h>
#include <vtkQuadraticTriangle.h>
#include <vtkQuadraticQuad.h>
#include <vtkTetra.h>
#include <vtkPyramid.h>
#include <vtkWedge.h>
#include <vtkHexahedron.h>
#include <vtkQuadraticTetra.h>
#include <vtkQuadraticPyramid.h>
#include <vtkQuadraticWedge.h>
#include <vtkQuadraticHexahedron.h>
#include "FemVTKTools.h"
#include "FemMeshProperty.h"
@@ -212,6 +214,7 @@ FemMesh* FemVTKTools::readVTKMesh(const char* filename, FemMesh* mesh)
void exportFemMeshFaces(vtkSmartPointer<vtkUnstructuredGrid> grid, const SMDS_FaceIteratorPtr& aFaceIter)
{
Base::Console().Log(" Start: VTK mesh builder faces.\n");
vtkSmartPointer<vtkCellArray> triangleArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> quadTriangleArray = vtkSmartPointer<vtkCellArray>::New();
@@ -274,57 +277,60 @@ void exportFemMeshFaces(vtkSmartPointer<vtkUnstructuredGrid> 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<vtkUnstructuredGrid> 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<vtkCellArray> tetraArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> hexaArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> wedgeArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> pyramidArray = vtkSmartPointer<vtkCellArray>::New();
// quadratic elemnts with 13 and 15 nodes are not added yet
vtkSmartPointer<vtkCellArray> wedgeArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> hexaArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> quadTetraArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> quadPyramidArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> quadWedgeArray = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> quadHexaArray = vtkSmartPointer<vtkCellArray>::New();
for (;aVolIter->more();)
{
const SMDS_MeshVolume* aVol = aVolIter->next();
if (aVol->NbNodes() == 4) { // tetrahedra
vtkSmartPointer<vtkTetra> tetra = vtkSmartPointer<vtkTetra>::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<vtkTetra> cell = vtkSmartPointer<vtkTetra>::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<vtkPyramid> cell= vtkSmartPointer<vtkPyramid>::New();
else if (aVol->NbNodes() == 5) { // pyra5
Base::Console().Log(" Volume pyra5\n");
vtkSmartPointer<vtkPyramid> cell = vtkSmartPointer<vtkPyramid>::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<vtkWedge> cell = vtkSmartPointer<vtkWedge>::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<vtkUnstructuredGrid> 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<vtkHexahedron> cell= vtkSmartPointer<vtkHexahedron>::New();
else if (aVol->NbNodes() == 8) { // hexa8
Base::Console().Log(" Volume hexa8\n");
vtkSmartPointer<vtkHexahedron> cell = vtkSmartPointer<vtkHexahedron>::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<vtkUnstructuredGrid> 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<vtkQuadraticTetra> tetra = vtkSmartPointer<vtkQuadraticTetra>::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<vtkHexahedron> cell= vtkSmartPointer<vtkHexahedron>::New();
else if (aVol->NbNodes() == 13) { // pyra13
Base::Console().Log(" Volume pyra13\n");
vtkSmartPointer<vtkQuadraticPyramid> cell = vtkSmartPointer<vtkQuadraticPyramid>::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<vtkQuadraticWedge> cell = vtkSmartPointer<vtkQuadraticWedge>::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<vtkQuadraticHexahedron> cell = vtkSmartPointer<vtkQuadraticHexahedron>::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<vtkUnstructuredGrid> 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<vtkUnstructuredGrid> grid, float scale)
{
Base::Console().Message("Start: VTK mesh builder ======================\n");
SMESH_Mesh* smesh = const_cast<SMESH_Mesh*>(mesh->getSMesh());
SMESHDS_Mesh* meshDS = smesh->GetMeshDS();
const SMDS_MeshInfo& info = meshDS->GetMeshInfo();
@@ -410,20 +438,21 @@ void FemVTKTools::exportVTKMesh(const FemMesh* mesh, vtkSmartPointer<vtkUnstruct
points->SetPoint(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<vtkUnstructuredGrid> grid = vtkSmartPointer<vtkUnstructuredGrid>::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<PropertyFemMesh*>(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