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

Fem: Export/import mesh elements to Vtk using proper order

Browse Source
This commit is contained in:
marioalexis
2024-03-09 15:05:53 -03:00
parent 254ed25e07
commit 1225e7be17
2 changed files with 157 additions and 237 deletions
Download Patch File Download Diff File Expand all files Collapse all files

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

@@ -97,6 +97,52 @@ void writeVTKFile(const char* filename, vtkSmartPointer<vtkUnstructuredGrid> dat
writer->Write();
}
namespace
{
// Helper function to fill vtkCellArray from SMDS_Mesh using vtk cell order
template<typename T, typename E>
void fillVtkArray(vtkSmartPointer<vtkCellArray>& elemArray, std::vector<int>& types, const E* elem)
{
vtkSmartPointer<T> cell = vtkSmartPointer<T>::New();
const std::vector<int>& order = SMDS_MeshCell::toVtkOrder(elem->GetEntityType());
if (!order.empty()) {
for (int i = 0; i < elem->NbNodes(); ++i) {
cell->GetPointIds()->SetId(i, elem->GetNode(order[i])->GetID() - 1);
}
}
else {
for (int i = 0; i < elem->NbNodes(); ++i) {
cell->GetPointIds()->SetId(i, elem->GetNode(i)->GetID() - 1);
}
}
elemArray->InsertNextCell(cell);
types.push_back(SMDS_MeshCell::toVtkType(elem->GetEntityType()));
}
// Helper function to fill SMDS_Mesh elements ID from vtk cell
void fillMeshElementIds(vtkCell* cell, std::vector<int>& ids)
{
VTKCellType cellType = static_cast<VTKCellType>(cell->GetCellType());
const std::vector<int>& order = SMDS_MeshCell::fromVtkOrder(cellType);
vtkIdType* vtkIds = cell->GetPointIds()->GetPointer(0);
ids.clear();
int nbPoints = cell->GetNumberOfPoints();
ids.resize(nbPoints);
if (!order.empty()) {
for (int i = 0; i < nbPoints; ++i) {
ids[i] = vtkIds[order[i]] + 1;
}
}
else {
for (int i = 0; i < nbPoints; ++i) {
ids[i] = vtkIds[i] + 1;
}
}
}
} // namespace
void FemVTKTools::importVTKMesh(vtkSmartPointer<vtkDataSet> dataset, FemMesh* mesh, float scale)
{
@@ -105,10 +151,6 @@ void FemVTKTools::importVTKMesh(vtkSmartPointer<vtkDataSet> dataset, FemMesh* me
Base::Console().Log("%d nodes/points and %d cells/elements found!\n", nPoints, nCells);
Base::Console().Log("Build SMESH mesh out of the vtk mesh data.\n", nPoints, nCells);
// vtkSmartPointer<vtkCellArray> cells = dataset->GetCells();
// works only for vtkUnstructuredGrid
vtkSmartPointer<vtkIdList> idlist = vtkSmartPointer<vtkIdList>::New();
// Now fill the SMESH datastructure
SMESH_Mesh* smesh = mesh->getSMesh();
SMESHDS_Mesh* meshds = smesh->GetMeshDS();
@@ -120,138 +162,120 @@ void FemVTKTools::importVTKMesh(vtkSmartPointer<vtkDataSet> dataset, FemMesh* me
}
for (vtkIdType iCell = 0; iCell < nCells; iCell++) {
idlist->Reset();
idlist = dataset->GetCell(iCell)->GetPointIds();
vtkIdType* ids = idlist->GetPointer(0);
switch (dataset->GetCellType(iCell)) {
vtkCell* cell = dataset->GetCell(iCell);
std::vector<int> ids;
fillMeshElementIds(cell, ids);
switch (cell->GetCellType()) {
// 2D faces
case VTK_TRIANGLE: // tria3
meshds->AddFaceWithID(ids[0] + 1, ids[1] + 1, ids[2] + 1, iCell + 1);
meshds->AddFaceWithID(ids[0], ids[1], ids[2], iCell + 1);
break;
case VTK_QUADRATIC_TRIANGLE: // tria6
meshds->AddFaceWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
iCell + 1);
meshds->AddFaceWithID(ids[0], ids[1], ids[2], ids[3], ids[4], ids[5], iCell + 1);
break;
case VTK_QUAD: // quad4
meshds->AddFaceWithID(ids[0] + 1, ids[1] + 1, ids[2] + 1, ids[3] + 1, iCell + 1);
meshds->AddFaceWithID(ids[0], ids[1], ids[2], ids[3], iCell + 1);
break;
case VTK_QUADRATIC_QUAD: // quad8
meshds->AddFaceWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
ids[6] + 1,
ids[7] + 1,
meshds->AddFaceWithID(ids[0],
ids[1],
ids[2],
ids[3],
ids[4],
ids[5],
ids[6],
ids[7],
iCell + 1);
break;
// 3D volumes
case VTK_TETRA: // tetra4
meshds->AddVolumeWithID(ids[0] + 1, ids[1] + 1, ids[2] + 1, ids[3] + 1, iCell + 1);
meshds->AddVolumeWithID(ids[0], ids[1], ids[2], ids[3], iCell + 1);
break;
case VTK_QUADRATIC_TETRA: // tetra10
meshds->AddVolumeWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
ids[6] + 1,
ids[7] + 1,
ids[8] + 1,
ids[9] + 1,
meshds->AddVolumeWithID(ids[0],
ids[1],
ids[2],
ids[3],
ids[4],
ids[5],
ids[6],
ids[7],
ids[8],
ids[9],
iCell + 1);
break;
case VTK_HEXAHEDRON: // hexa8
meshds->AddVolumeWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
ids[6] + 1,
ids[7] + 1,
meshds->AddVolumeWithID(ids[0],
ids[1],
ids[2],
ids[3],
ids[4],
ids[5],
ids[6],
ids[7],
iCell + 1);
break;
case VTK_QUADRATIC_HEXAHEDRON: // hexa20
meshds->AddVolumeWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
ids[6] + 1,
ids[7] + 1,
ids[8] + 1,
ids[9] + 1,
ids[10] + 1,
ids[11] + 1,
ids[12] + 1,
ids[13] + 1,
ids[14] + 1,
ids[15] + 1,
ids[16] + 1,
ids[17] + 1,
ids[18] + 1,
ids[19] + 1,
meshds->AddVolumeWithID(ids[0],
ids[1],
ids[2],
ids[3],
ids[4],
ids[5],
ids[6],
ids[7],
ids[8],
ids[9],
ids[10],
ids[11],
ids[12],
ids[13],
ids[14],
ids[15],
ids[16],
ids[17],
ids[18],
ids[19],
iCell + 1);
break;
case VTK_WEDGE: // penta6
meshds->AddVolumeWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
iCell + 1);
meshds->AddVolumeWithID(ids[0], ids[1], ids[2], ids[3], ids[4], ids[5], iCell + 1);
break;
case VTK_QUADRATIC_WEDGE: // penta15
meshds->AddVolumeWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
ids[6] + 1,
ids[7] + 1,
ids[8] + 1,
ids[9] + 1,
ids[10] + 1,
ids[11] + 1,
ids[12] + 1,
ids[13] + 1,
ids[14] + 1,
meshds->AddVolumeWithID(ids[0],
ids[1],
ids[2],
ids[3],
ids[4],
ids[5],
ids[6],
ids[7],
ids[8],
ids[9],
ids[10],
ids[11],
ids[12],
ids[13],
ids[14],
iCell + 1);
break;
case VTK_PYRAMID: // pyra5
meshds->AddVolumeWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
iCell + 1);
meshds->AddVolumeWithID(ids[0], ids[1], ids[2], ids[3], ids[4], iCell + 1);
break;
case VTK_QUADRATIC_PYRAMID: // pyra13
meshds->AddVolumeWithID(ids[0] + 1,
ids[1] + 1,
ids[2] + 1,
ids[3] + 1,
ids[4] + 1,
ids[5] + 1,
ids[6] + 1,
ids[7] + 1,
ids[8] + 1,
ids[9] + 1,
ids[10] + 1,
ids[11] + 1,
ids[12] + 1,
meshds->AddVolumeWithID(ids[0],
ids[1],
ids[2],
ids[3],
ids[4],
ids[5],
ids[6],
ids[7],
ids[8],
ids[9],
ids[10],
ids[11],
ids[12],
iCell + 1);
break;
@@ -314,58 +338,23 @@ void exportFemMeshFaces(vtkSmartPointer<vtkUnstructuredGrid> grid,
vtkSmartPointer<vtkCellArray> elemArray = vtkSmartPointer<vtkCellArray>::New();
std::vector<int> types;
for (; aFaceIter->more();) {
while (aFaceIter->more()) {
const SMDS_MeshFace* aFace = aFaceIter->next();
// triangle
if (aFace->NbNodes() == 3) {
vtkSmartPointer<vtkTriangle> tria = vtkSmartPointer<vtkTriangle>::New();
tria->GetPointIds()->SetId(0, aFace->GetNode(0)->GetID() - 1);
tria->GetPointIds()->SetId(1, aFace->GetNode(1)->GetID() - 1);
tria->GetPointIds()->SetId(2, aFace->GetNode(2)->GetID() - 1);
elemArray->InsertNextCell(tria);
types.push_back(VTK_TRIANGLE);
if (aFace->GetEntityType() == SMDSEntity_Triangle) {
fillVtkArray<vtkTriangle>(elemArray, types, aFace);
}
// quad
else if (aFace->NbNodes() == 4) {
vtkSmartPointer<vtkQuad> quad = vtkSmartPointer<vtkQuad>::New();
quad->GetPointIds()->SetId(0, aFace->GetNode(0)->GetID() - 1);
quad->GetPointIds()->SetId(1, aFace->GetNode(1)->GetID() - 1);
quad->GetPointIds()->SetId(2, aFace->GetNode(2)->GetID() - 1);
quad->GetPointIds()->SetId(3, aFace->GetNode(3)->GetID() - 1);
elemArray->InsertNextCell(quad);
types.push_back(VTK_QUAD);
else if (aFace->GetEntityType() == SMDSEntity_Quadrangle) {
fillVtkArray<vtkQuad>(elemArray, types, aFace);
}
// quadratic triangle
else if (aFace->NbNodes() == 6) {
vtkSmartPointer<vtkQuadraticTriangle> tria =
vtkSmartPointer<vtkQuadraticTriangle>::New();
tria->GetPointIds()->SetId(0, aFace->GetNode(0)->GetID() - 1);
tria->GetPointIds()->SetId(1, aFace->GetNode(1)->GetID() - 1);
tria->GetPointIds()->SetId(2, aFace->GetNode(2)->GetID() - 1);
tria->GetPointIds()->SetId(3, aFace->GetNode(3)->GetID() - 1);
tria->GetPointIds()->SetId(4, aFace->GetNode(4)->GetID() - 1);
tria->GetPointIds()->SetId(5, aFace->GetNode(5)->GetID() - 1);
elemArray->InsertNextCell(tria);
types.push_back(VTK_QUADRATIC_TRIANGLE);
else if (aFace->GetEntityType() == SMDSEntity_Quad_Triangle) {
fillVtkArray<vtkQuadraticTriangle>(elemArray, types, aFace);
}
// quadratic quad
else if (aFace->NbNodes() == 8) {
vtkSmartPointer<vtkQuadraticQuad> quad = vtkSmartPointer<vtkQuadraticQuad>::New();
quad->GetPointIds()->SetId(0, aFace->GetNode(0)->GetID() - 1);
quad->GetPointIds()->SetId(1, aFace->GetNode(1)->GetID() - 1);
quad->GetPointIds()->SetId(2, aFace->GetNode(2)->GetID() - 1);
quad->GetPointIds()->SetId(3, aFace->GetNode(3)->GetID() - 1);
quad->GetPointIds()->SetId(4, aFace->GetNode(4)->GetID() - 1);
quad->GetPointIds()->SetId(5, aFace->GetNode(5)->GetID() - 1);
quad->GetPointIds()->SetId(6, aFace->GetNode(6)->GetID() - 1);
quad->GetPointIds()->SetId(7, aFace->GetNode(7)->GetID() - 1);
elemArray->InsertNextCell(quad);
types.push_back(VTK_QUADRATIC_QUAD);
else if (aFace->GetEntityType() == SMDSEntity_Quad_Quadrangle) {
fillVtkArray<vtkQuadraticQuad>(elemArray, types, aFace);
}
else {
throw Base::TypeError("Face not yet supported by FreeCAD's VTK mesh builder\n");
@@ -387,101 +376,32 @@ void exportFemMeshCells(vtkSmartPointer<vtkUnstructuredGrid> grid,
vtkSmartPointer<vtkCellArray> elemArray = vtkSmartPointer<vtkCellArray>::New();
std::vector<int> types;
for (; aVolIter->more();) {
while (aVolIter->more()) {
const SMDS_MeshVolume* aVol = aVolIter->next();
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);
elemArray->InsertNextCell(cell);
types.push_back(VTK_TETRA);
if (aVol->GetEntityType() == SMDSEntity_Tetra) { // tetra4
fillVtkArray<vtkTetra>(elemArray, types, aVol);
}
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);
elemArray->InsertNextCell(cell);
types.push_back(VTK_PYRAMID);
else if (aVol->GetEntityType() == SMDSEntity_Pyramid) { // pyra5
fillVtkArray<vtkPyramid>(elemArray, types, aVol);
}
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);
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);
cell->GetPointIds()->SetId(5, aVol->GetNode(5)->GetID() - 1);
elemArray->InsertNextCell(cell);
types.push_back(VTK_WEDGE);
else if (aVol->GetEntityType() == SMDSEntity_Penta) { // penta6
fillVtkArray<vtkWedge>(elemArray, types, aVol);
}
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);
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);
cell->GetPointIds()->SetId(6, aVol->GetNode(6)->GetID() - 1);
cell->GetPointIds()->SetId(7, aVol->GetNode(7)->GetID() - 1);
elemArray->InsertNextCell(cell);
types.push_back(VTK_HEXAHEDRON);
else if (aVol->GetEntityType() == SMDSEntity_Hexa) { // hexa8
fillVtkArray<vtkHexahedron>(elemArray, types, aVol);
}
else if (aVol->NbNodes() == 10) { // tetra10
Base::Console().Log(" Volume tetra10\n");
vtkSmartPointer<vtkQuadraticTetra> cell = vtkSmartPointer<vtkQuadraticTetra>::New();
for (int i = 0; i < 10; i++) {
cell->GetPointIds()->SetId(i, aVol->GetNode(i)->GetID() - 1);
}
elemArray->InsertNextCell(cell);
types.push_back(VTK_QUADRATIC_TETRA);
else if (aVol->GetEntityType() == SMDSEntity_Quad_Tetra) { // tetra10
fillVtkArray<vtkQuadraticTetra>(elemArray, types, aVol);
}
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);
}
elemArray->InsertNextCell(cell);
types.push_back(VTK_QUADRATIC_PYRAMID);
else if (aVol->GetEntityType() == SMDSEntity_Quad_Pyramid) { // pyra13
fillVtkArray<vtkQuadraticPyramid>(elemArray, types, aVol);
}
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);
}
elemArray->InsertNextCell(cell);
types.push_back(VTK_QUADRATIC_WEDGE);
else if (aVol->GetEntityType() == SMDSEntity_Quad_Penta) { // penta15
fillVtkArray<vtkQuadraticWedge>(elemArray, types, aVol);
}
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);
}
elemArray->InsertNextCell(cell);
types.push_back(VTK_QUADRATIC_HEXAHEDRON);
else if (aVol->GetEntityType() == SMDSEntity_Quad_Hexa) { // hexa20
fillVtkArray<vtkQuadraticHexahedron>(elemArray, types, aVol);
}
else {
throw Base::TypeError("Volume not yet supported by FreeCAD's VTK mesh builder\n");