1 // Copyright (C) 2010-2021 CEA/DEN, EDF R&D, OPEN CASCADE
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 // File: SMDS_Downward.cxx
21 // Created: Jun 3, 2010
24 #include "SMDS_Downward.hxx"
25 #include "SMDS_Mesh.hxx"
26 #include "utilities.h"
28 #include <vtkCellType.h>
29 #include <vtkCellLinks.h>
30 #include <smIdType.hxx>
36 // ---------------------------------------------------------------------------
38 vector<int> SMDS_Downward::_cellDimension;
40 /*! get the dimension of a cell (1,2,3 for 1D, 2D 3D) given the vtk cell type
42 * @param cellType vtk cell type @see vtkCellType.h
45 int SMDS_Downward::getCellDimension(unsigned char cellType)
47 if (_cellDimension.empty())
49 _cellDimension.resize(VTK_MAXTYPE + 1, 0);
50 _cellDimension[VTK_LINE] = 1;
51 _cellDimension[VTK_QUADRATIC_EDGE] = 1;
52 _cellDimension[VTK_TRIANGLE] = 2;
53 _cellDimension[VTK_QUADRATIC_TRIANGLE] = 2;
54 _cellDimension[VTK_BIQUADRATIC_TRIANGLE] = 2;
55 _cellDimension[VTK_QUAD] = 2;
56 _cellDimension[VTK_QUADRATIC_QUAD] = 2;
57 _cellDimension[VTK_BIQUADRATIC_QUAD] = 2;
58 _cellDimension[VTK_TETRA] = 3;
59 _cellDimension[VTK_QUADRATIC_TETRA] = 3;
60 _cellDimension[VTK_HEXAHEDRON] = 3;
61 _cellDimension[VTK_QUADRATIC_HEXAHEDRON] = 3;
62 _cellDimension[VTK_TRIQUADRATIC_HEXAHEDRON] = 3;
63 _cellDimension[VTK_WEDGE] = 3;
64 _cellDimension[VTK_QUADRATIC_WEDGE] = 3;
65 _cellDimension[VTK_BIQUADRATIC_QUADRATIC_WEDGE] = 3;
66 _cellDimension[VTK_PYRAMID] = 3;
67 _cellDimension[VTK_QUADRATIC_PYRAMID] = 3;
68 _cellDimension[VTK_HEXAGONAL_PRISM] = 3;
70 return _cellDimension[cellType];
73 // ---------------------------------------------------------------------------
75 /*! Generic constructor for all the downward connectivity structures (one per vtk cell type).
76 * The static structure for cell dimension is set only once.
77 * @param grid unstructured grid associated to the mesh.
78 * @param nbDownCells number of downward entities associated to this vtk type of cell.
81 SMDS_Downward::SMDS_Downward(SMDS_UnstructuredGrid *grid, int nbDownCells) :
82 _grid(grid), _nbDownCells(nbDownCells)
85 this->_cellIds.clear();
86 this->_cellTypes.clear();
87 if (_cellDimension.empty())
88 getCellDimension( VTK_LINE );
91 SMDS_Downward::~SMDS_Downward()
95 /*! Give or create an entry for downward connectivity structure relative to a cell.
96 * If the entry already exists, just return its id, otherwise, create it.
97 * The internal storage memory is allocated if needed.
98 * The SMDS_UnstructuredGrid::_cellIdToDownId vector is completed for vtkUnstructuredGrid cells.
99 * @param vtkId for a vtkUnstructuredGrid cell or -1 (default) for a created downward cell.
100 * @return the rank in downward[vtkType] structure.
102 int SMDS_Downward::addCell(int vtkId)
106 localId = _grid->CellIdToDownId(vtkId);
110 localId = this->_maxId;
112 this->allocate(_maxId);
115 this->_vtkCellIds[localId] = vtkId;
116 _grid->setCellIdToDownId(vtkId, localId);
118 this->initCell(localId);
122 /*! generic method do nothing. see derived methods
126 void SMDS_Downward::initCell(int /*cellId*/)
130 /*! Get the number of downward entities associated to a cell (always the same for a given vtk type of cell)
132 * @param cellId not used here.
135 int SMDS_Downward::getNumberOfDownCells(int /*cellId*/)
140 /*! get a pointer on the downward entities id's associated to a cell.
141 * @see SMDS_Downward::getNumberOfDownCells for the number of downward entities.
142 * @see SMDS_Downward::getDownTypes for the vtk cell types associated to the downward entities.
143 * @param cellId index of the cell in the downward structure relative to a given vtk cell type.
144 * @return table of downward entities id's.
146 const int* SMDS_Downward::getDownCells(int cellId)
148 //ASSERT((cellId >=0) && (cellId < _maxId));
149 return &_cellIds[_nbDownCells * cellId];
152 /*! get a list of vtk cell types associated to downward entities of a given cell, in the same order
153 * than the downward entities id's list (@see SMDS_Downward::getDownCells).
155 * @param cellId index of the cell in the downward structure relative to a vtk cell type.
156 * @return table of downward entities types.
158 const unsigned char* SMDS_Downward::getDownTypes(int /*cellId*/)
160 return &_cellTypes[0];
163 /*! add a downward entity of dimension n-1 (cell or node) to a given cell.
164 * Actual implementation is done in derived methods.
165 * @param cellId index of the parent cell (dimension n) in the downward structure relative to a vtk cell type.
166 * @param lowCellId index of the children cell to add (dimension n-1)
167 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
169 void SMDS_Downward::addDownCell(int /*cellId*/, int /*lowCellId*/, unsigned char /*aType*/)
171 ASSERT(0); // must be re-implemented in derived class
174 /*! add a downward entity of dimension n+1 to a given cell.
175 * Actual implementation is done in derived methods.
176 * @param cellId index of the children cell (dimension n) in the downward structure relative to a vtk cell type.
177 * @param upCellId index of the parent cell to add (dimension n+1)
178 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
180 void SMDS_Downward::addUpCell(int /*cellId*/, int /*upCellId*/, unsigned char /*aType*/)
182 ASSERT(0); // must be re-implemented in derived class
185 int SMDS_Downward::getNodeSet(int /*cellId*/, int* /*nodeSet*/)
190 // ---------------------------------------------------------------------------
192 SMDS_Down1D::SMDS_Down1D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
193 SMDS_Downward(grid, nbDownCells)
195 _upCellIdsVector.clear();
196 _upCellTypesVector.clear();
198 _upCellTypes.clear();
199 _upCellIndex.clear();
202 SMDS_Down1D::~SMDS_Down1D()
206 /*! clear vectors used to reference 2D cells containing the edge
210 void SMDS_Down1D::initCell(int cellId)
212 _upCellIdsVector[cellId].clear();
213 _upCellTypesVector[cellId].clear();
216 /*! Resize the downward connectivity storage vector if needed.
218 * @param nbElems total number of elements of the same type required
220 void SMDS_Down1D::allocate(int nbElems)
222 if (nbElems >= (int)_vtkCellIds.size())
224 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
225 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
226 _upCellIdsVector.resize(nbElems + SMDS_Mesh::chunkSize);
227 _upCellTypesVector.resize(nbElems + SMDS_Mesh::chunkSize);
231 void SMDS_Down1D::compactStorage()
233 _cellIds.resize(_nbDownCells * _maxId);
234 _vtkCellIds.resize(_maxId);
236 smIdType sizeUpCells = 0;
237 for (int i = 0; i < _maxId; i++)
238 sizeUpCells += _upCellIdsVector[i].size();
239 _upCellIds.resize(sizeUpCells, -1);
240 _upCellTypes.resize(sizeUpCells);
241 _upCellIndex.resize(_maxId + 1, -1); // id and types of rank i correspond to [ _upCellIndex[i], _upCellIndex[i+1] [
243 for (int i = 0; i < _maxId; i++)
245 _upCellIndex[i] = current;
246 for (size_t j = 0; j < _upCellIdsVector[i].size(); j++)
248 _upCellIds[current] = _upCellIdsVector[i][j];
249 _upCellTypes[current] = _upCellTypesVector[i][j];
253 _upCellIndex[_maxId] = current;
255 _upCellIdsVector.clear();
256 _upCellTypesVector.clear();
259 void SMDS_Down1D::addUpCell(int cellId, int upCellId, unsigned char aType)
261 //ASSERT((cellId >=0) && (cellId < _maxId));
262 smIdType nbFaces = _upCellIdsVector[cellId].size();
263 for (smIdType i = 0; i < nbFaces; i++)
265 if ((_upCellIdsVector[cellId][i] == upCellId) && (_upCellTypesVector[cellId][i] == aType))
267 return; // already done
270 _upCellIdsVector[cellId].push_back(upCellId);
271 _upCellTypesVector[cellId].push_back(aType);
274 int SMDS_Down1D::getNumberOfUpCells(int cellId)
276 //ASSERT((cellId >=0) && (cellId < _maxId));
277 return _upCellIndex[cellId + 1] - _upCellIndex[cellId];
280 const int* SMDS_Down1D::getUpCells(int cellId)
282 //ASSERT((cellId >=0) && (cellId < _maxId));
283 return &_upCellIds[_upCellIndex[cellId]];
286 const unsigned char* SMDS_Down1D::getUpTypes(int cellId)
288 //ASSERT((cellId >=0) && (cellId < _maxId));
289 return &_upCellTypes[_upCellIndex[cellId]];
292 void SMDS_Down1D::getNodeIds(int cellId, std::set<int>& nodeSet)
294 for (int i = 0; i < _nbDownCells; i++)
295 nodeSet.insert(_cellIds[_nbDownCells * cellId + i]);
298 int SMDS_Down1D::getNodeSet(int cellId, int* nodeSet)
300 for (int i = 0; i < _nbDownCells; i++)
301 nodeSet[i] = _cellIds[_nbDownCells * cellId + i];
305 void SMDS_Down1D::setNodes(int cellId, int vtkId)
308 vtkIdType const *pts(nullptr); // will refer to the point id's of the face
309 _grid->GetCellPoints(vtkId, npts, pts);
310 // MESSAGE(vtkId << " " << npts << " " << _nbDownCells);
311 //ASSERT(npts == _nbDownCells);
312 for (int i = 0; i < npts; i++)
314 _cellIds[_nbDownCells * cellId + i] = pts[i];
318 void SMDS_Down1D::setNodes(int cellId, const int* nodeIds)
320 //ASSERT(nodeIds.size() == _nbDownCells);
321 for (int i = 0; i < _nbDownCells; i++)
323 _cellIds[_nbDownCells * cellId + i] = nodeIds[i];
327 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
328 * We keep in the list the cells that contains all the nodes, we keep only volumes and faces.
329 * @param cellId id of the edge in the downward structure
330 * @param vtkIds vector of vtk id's
331 * @return number of vtk cells (size of vector)
333 int SMDS_Down1D::computeVtkCells(int cellId, std::vector<int>& vtkIds)
337 // --- find all the cells the points belong to, and how many of the points belong to a given cell
339 int *pts = &_cellIds[_nbDownCells * cellId];
340 int ncells = this->computeVtkCells(pts, vtkIds);
344 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
346 * @param pts list of points id's defining an edge
347 * @param vtkIds vector of vtk id's
348 * @return number of vtk cells (size of vector)
350 int SMDS_Down1D::computeVtkCells(int *pts, std::vector<int>& vtkIds)
353 // --- find all the cells the points belong to, and how many of the points belong to a given cell
358 for (int i = 0; i < _nbDownCells; i++)
360 vtkIdType point = pts[i];
361 int numCells = _grid->GetLinks()->GetNcells(point);
362 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
363 for (int j = 0; j < numCells; j++)
365 int vtkCellId = cells[j];
367 for (int k = 0; k < cnt; k++)
369 if (cellIds[k] == vtkCellId)
378 cellIds[cnt] = vtkCellId;
380 // TODO ASSERT(cnt<1000);
386 // --- find the face and volume cells: they contains all the points and are of type volume or face
389 for (int i = 0; i < cnt; i++)
391 if (cellCnt[i] == _nbDownCells)
393 int vtkElemId = cellIds[i];
394 int vtkType = _grid->GetCellType(vtkElemId);
395 if (SMDS_Downward::getCellDimension(vtkType) > 1)
397 vtkIds.push_back(vtkElemId);
406 /*! Build the list of downward faces from a list of vtk cells.
408 * @param cellId id of the edge in the downward structure
409 * @param vtkIds vector of vtk id's
410 * @param downFaces vector of face id's in downward structures
411 * @param downTypes vector of face types
412 * @return number of downward faces
414 int SMDS_Down1D::computeFaces(int cellId, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
416 int *pts = &_cellIds[_nbDownCells * cellId];
417 int nbFaces = this->computeFaces(pts, vtkIds, nbcells, downFaces, downTypes);
421 /*! Build the list of downward faces from a list of vtk cells.
423 * @param pts list of points id's defining an edge
424 * @param vtkIds vector of vtk id's
425 * @param downFaces vector of face id's in downward structures
426 * @param downTypes vector of face types
427 * @return number of downward faces
429 int SMDS_Down1D::computeFaces(int* pts, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
432 for (int i = 0; i < nbcells; i++)
434 int vtkId = vtkIds[i];
435 int vtkType = _grid->GetCellType(vtkId);
436 if (SMDS_Downward::getCellDimension(vtkType) == 2)
438 int faceId = _grid->CellIdToDownId(vtkId);
439 downFaces[cnt] = faceId;
440 downTypes[cnt] = vtkType;
443 else if (SMDS_Downward::getCellDimension(vtkType) == 3)
445 int volId = _grid->CellIdToDownId(vtkId);
446 SMDS_Downward * downvol = _grid->getDownArray(vtkType);
447 //const int *downIds = downvol->getDownCells(volId);
448 const unsigned char* downTypesVol = downvol->getDownTypes(volId);
449 int nbFaces = downvol->getNumberOfDownCells(volId);
450 const int* faceIds = downvol->getDownCells(volId);
451 for (int n = 0; n < nbFaces; n++)
453 SMDS_Down2D *downFace = static_cast<SMDS_Down2D*> (_grid->getDownArray(downTypesVol[n]));
454 bool isInFace = downFace->isInFace(faceIds[n], pts, _nbDownCells);
457 bool alreadySet = false;
458 for (int k = 0; k < cnt; k++)
459 if (faceIds[n] == downFaces[k])
466 downFaces[cnt] = faceIds[n];
467 downTypes[cnt] = downTypesVol[n];
477 // ---------------------------------------------------------------------------
479 SMDS_Down2D::SMDS_Down2D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
480 SMDS_Downward(grid, nbDownCells)
483 _upCellTypes.clear();
488 SMDS_Down2D::~SMDS_Down2D()
492 int SMDS_Down2D::getNumberOfUpCells(int cellId)
495 if (_upCellIds[2 * cellId] >= 0)
497 if (_upCellIds[2 * cellId + 1] >= 0)
502 const int* SMDS_Down2D::getUpCells(int cellId)
504 //ASSERT((cellId >=0) && (cellId < _maxId));
505 return &_upCellIds[2 * cellId];
508 const unsigned char* SMDS_Down2D::getUpTypes(int cellId)
510 //ASSERT((cellId >=0) && (cellId < _maxId));
511 return &_upCellTypes[2 * cellId];
514 void SMDS_Down2D::getNodeIds(int cellId, std::set<int>& nodeSet)
516 for (int i = 0; i < _nbDownCells; i++)
518 int downCellId = _cellIds[_nbDownCells * cellId + i];
519 unsigned char cellType = _cellTypes[i];
520 this->_grid->getDownArray(cellType)->getNodeIds(downCellId, nodeSet);
524 /*! Find in vtkUnstructuredGrid the volumes containing a face already stored in vtkUnstructuredGrid.
525 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
526 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
527 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
528 * @param cellId the face cell id in vkUnstructuredGrid
529 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
530 * @return number of volumes (0, 1 or 2)
532 int SMDS_Down2D::computeVolumeIds(int cellId, int* ids)
534 // --- find point id's of the face
537 vtkIdType const *pts(nullptr); // will refer to the point id's of the face
538 _grid->GetCellPoints(cellId, npts, pts);
540 for (int i = 0; i < npts; i++)
541 nodes.push_back(pts[i]);
542 int nvol = this->computeVolumeIdsFromNodesFace(&nodes[0], npts, ids);
546 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
547 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
548 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
549 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
551 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
552 * @return number of volumes (0, 1 or 2)
554 int SMDS_Down2D::computeVolumeIds(ElemByNodesType& faceByNodes, int* ids)
556 int nvol = this->computeVolumeIdsFromNodesFace(&faceByNodes.nodeIds[0], faceByNodes.nbNodes, ids);
560 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
561 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
562 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
563 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
564 * @param pts array of vtk node id's
565 * @param npts number of nodes
567 * @return number of volumes (0, 1 or 2)
569 int SMDS_Down2D::computeVolumeIdsFromNodesFace(int* pts, int npts, int* ids)
572 // --- find all the cells the points belong to, and how many of the points belong to a given cell
577 for (int i = 0; i < npts; i++)
579 vtkIdType point = pts[i];
580 int numCells = _grid->GetLinks()->GetNcells(point);
581 //MESSAGE("cells pour " << i << " " << numCells);
582 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
583 for (int j = 0; j < numCells; j++)
585 int vtkCellId = cells[j];
587 for (int k = 0; k < cnt; k++)
589 if (cellIds[k] == vtkCellId)
598 cellIds[cnt] = vtkCellId;
600 // TODO ASSERT(cnt<1000);
606 // --- find the volume cells: they contains all the points and are of type volume
609 for (int i = 0; i < cnt; i++)
611 //MESSAGE("cell " << cellIds[i] << " points " << cellCnt[i]);
612 if (cellCnt[i] == npts)
614 int vtkElemId = cellIds[i];
615 int vtkType = _grid->GetCellType(vtkElemId);
616 if (SMDS_Downward::getCellDimension(vtkType) == 3)
618 ids[nvol] = vtkElemId; // store the volume id in given vector
629 void SMDS_Down2D::setTempNodes(int cellId, int vtkId)
632 vtkIdType const *pts(nullptr); // will refer to the point id's of the face
633 _grid->GetCellPoints(vtkId, npts, pts);
634 // MESSAGE(vtkId << " " << npts << " " << _nbNodes);
635 //ASSERT(npts == _nbNodes);
636 for (int i = 0; i < npts; i++)
638 _tempNodes[_nbNodes * cellId + i] = pts[i];
642 void SMDS_Down2D::setTempNodes(int cellId, ElemByNodesType& faceByNodes)
644 for (int i = 0; i < faceByNodes.nbNodes; i++)
645 _tempNodes[_nbNodes * cellId + i] = faceByNodes.nodeIds[i];
648 /*! Find if all the nodes belongs to the face.
650 * @param cellId the face cell Id
651 * @param nodeSet set of node id's to be found in the face list of nodes
654 bool SMDS_Down2D::isInFace(int cellId, int *pts, int npts)
657 int *nodes = &_tempNodes[_nbNodes * cellId];
658 for (int j = 0; j < npts; j++)
661 for (int i = 0; i < _nbNodes; i++)
663 if (nodes[i] == point)
670 return (nbFound == npts);
673 /*! Resize the downward connectivity storage vector if needed.
675 * @param nbElems total number of elements of the same type required
677 void SMDS_Down2D::allocate(int nbElems)
679 if (nbElems >= (int)_vtkCellIds.size())
681 _cellIds.resize (_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
682 _vtkCellIds.resize (nbElems + SMDS_Mesh::chunkSize, -1);
683 _upCellIds.resize (2 * (nbElems + SMDS_Mesh::chunkSize), -1);
684 _upCellTypes.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
685 _tempNodes.resize (_nbNodes * (nbElems + SMDS_Mesh::chunkSize), -1);
689 void SMDS_Down2D::compactStorage()
691 _cellIds.resize(_nbDownCells * _maxId);
692 _upCellIds.resize(2 * _maxId);
693 _upCellTypes.resize(2 * _maxId);
694 _vtkCellIds.resize(_maxId);
698 void SMDS_Down2D::addUpCell(int cellId, int upCellId, unsigned char aType)
700 //ASSERT((cellId >=0)&& (cellId < _maxId));
701 int *vols = &_upCellIds[2 * cellId];
702 unsigned char *types = &_upCellTypes[2 * cellId];
703 for (int i = 0; i < 2; i++)
707 vols[i] = upCellId; // use non affected volume
711 if ((vols[i] == upCellId) && (types[i] == aType)) // already done
717 int SMDS_Down2D::getNodeSet(int cellId, int* nodeSet)
719 for (int i = 0; i < _nbNodes; i++)
720 nodeSet[i] = _tempNodes[_nbNodes * cellId + i];
724 int SMDS_Down2D::FindEdgeByNodes(int cellId, ElemByNodesType& edgeByNodes)
726 int *edges = &_cellIds[_nbDownCells * cellId];
727 for (int i = 0; i < _nbDownCells; i++)
729 if ((edges[i] >= 0) && (edgeByNodes.vtkType == _cellTypes[i]))
732 int npts = this->_grid->getDownArray(edgeByNodes.vtkType)->getNodeSet(edges[i], nodeSet);
734 for (int j = 0; j < npts; j++)
736 int point = edgeByNodes.nodeIds[j];
738 for (int k = 0; k < npts; k++)
740 if (nodeSet[k] == point)
756 // ---------------------------------------------------------------------------
758 SMDS_Down3D::SMDS_Down3D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
759 SMDS_Downward(grid, nbDownCells)
763 SMDS_Down3D::~SMDS_Down3D()
767 void SMDS_Down3D::allocate(int nbElems)
769 if (nbElems >= (int)_vtkCellIds.size())
771 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
772 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
776 void SMDS_Down3D::compactStorage()
778 // nothing to do, size was known before
781 int SMDS_Down3D::getNumberOfUpCells(int /*cellId*/)
786 const int* SMDS_Down3D::getUpCells(int /*cellId*/)
791 const unsigned char* SMDS_Down3D::getUpTypes(int /*cellId*/)
796 void SMDS_Down3D::getNodeIds(int cellId, std::set<int>& nodeSet)
798 int vtkId = this->_vtkCellIds[cellId];
800 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
801 _grid->GetCellPoints(vtkId, npts, nodes);
802 for (int i = 0; i < npts; i++)
803 nodeSet.insert(nodes[i]);
806 int SMDS_Down3D::FindFaceByNodes(int cellId, ElemByNodesType& faceByNodes)
808 int *faces = &_cellIds[_nbDownCells * cellId];
811 for (int i = 0; i < _nbDownCells; i++)
813 if ((faces[i] >= 0) && (faceByNodes.vtkType == _cellTypes[i]))
816 npoints = faceByNodes.nbNodes;
819 int npts = this->_grid->getDownArray(faceByNodes.vtkType)->getNodeSet(faces[i], nodeSet);
821 continue; // skip this face
823 for (int j = 0; j < npts; j++)
825 int point = faceByNodes.nodeIds[j];
827 for (int k = 0; k < npts; k++)
829 if (nodeSet[k] == point)
832 break; // point j is in the 2 faces, skip remaining k values
836 break; // point j is not in the 2 faces, skip the remaining tests
845 // ---------------------------------------------------------------------------
847 SMDS_DownEdge::SMDS_DownEdge(SMDS_UnstructuredGrid *grid) :
850 _cellTypes.push_back(VTK_VERTEX);
851 _cellTypes.push_back(VTK_VERTEX);
854 SMDS_DownEdge::~SMDS_DownEdge()
858 // ---------------------------------------------------------------------------
860 SMDS_DownQuadEdge::SMDS_DownQuadEdge(SMDS_UnstructuredGrid *grid) :
863 _cellTypes.push_back(VTK_VERTEX);
864 _cellTypes.push_back(VTK_VERTEX);
865 _cellTypes.push_back(VTK_VERTEX);
868 SMDS_DownQuadEdge::~SMDS_DownQuadEdge()
872 // ---------------------------------------------------------------------------
874 SMDS_DownTriangle::SMDS_DownTriangle(SMDS_UnstructuredGrid *grid) :
877 _cellTypes.push_back(VTK_LINE);
878 _cellTypes.push_back(VTK_LINE);
879 _cellTypes.push_back(VTK_LINE);
883 SMDS_DownTriangle::~SMDS_DownTriangle()
887 void SMDS_DownTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
889 int *nodes = &_tempNodes[_nbNodes * cellId];
890 edgesWithNodes.nbElems = 3;
892 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
893 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
894 edgesWithNodes.elems[0].nbNodes = 2;
895 edgesWithNodes.elems[0].vtkType = VTK_LINE;
897 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
898 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
899 edgesWithNodes.elems[1].nbNodes = 2;
900 edgesWithNodes.elems[1].vtkType = VTK_LINE;
902 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
903 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
904 edgesWithNodes.elems[2].nbNodes = 2;
905 edgesWithNodes.elems[2].vtkType = VTK_LINE;
908 void SMDS_DownTriangle::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
910 //ASSERT((cellId >=0)&& (cellId < _maxId));
911 //ASSERT(aType == VTK_LINE);
912 int *faces = &_cellIds[_nbDownCells * cellId];
913 for (int i = 0; i < _nbDownCells; i++)
917 faces[i] = lowCellId;
920 if (faces[i] == lowCellId)
926 // ---------------------------------------------------------------------------
928 SMDS_DownQuadTriangle::SMDS_DownQuadTriangle(SMDS_UnstructuredGrid *grid) :
931 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
932 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
933 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
937 SMDS_DownQuadTriangle::~SMDS_DownQuadTriangle()
941 void SMDS_DownQuadTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
943 int *nodes = &_tempNodes[_nbNodes * cellId];
944 edgesWithNodes.nbElems = 3;
946 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
947 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
948 edgesWithNodes.elems[0].nodeIds[2] = nodes[3];
949 edgesWithNodes.elems[0].nbNodes = 3;
950 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
952 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
953 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
954 edgesWithNodes.elems[1].nodeIds[2] = nodes[4];
955 edgesWithNodes.elems[1].nbNodes = 3;
956 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
958 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
959 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
960 edgesWithNodes.elems[2].nodeIds[2] = nodes[5];
961 edgesWithNodes.elems[2].nbNodes = 3;
962 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
965 void SMDS_DownQuadTriangle::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
967 //ASSERT((cellId >=0)&& (cellId < _maxId));
968 //ASSERT(aType == VTK_QUADRATIC_EDGE);
969 int *edges = &_cellIds[_nbDownCells * cellId];
970 for (int i = 0; i < _nbDownCells; i++)
974 edges[i] = lowCellId;
977 if (edges[i] == lowCellId)
983 // ---------------------------------------------------------------------------
985 SMDS_DownQuadrangle::SMDS_DownQuadrangle(SMDS_UnstructuredGrid *grid) :
988 _cellTypes.push_back(VTK_LINE);
989 _cellTypes.push_back(VTK_LINE);
990 _cellTypes.push_back(VTK_LINE);
991 _cellTypes.push_back(VTK_LINE);
995 SMDS_DownQuadrangle::~SMDS_DownQuadrangle()
999 void SMDS_DownQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
1001 int *nodes = &_tempNodes[_nbNodes * cellId];
1002 edgesWithNodes.nbElems = 4;
1004 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1005 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1006 edgesWithNodes.elems[0].nbNodes = 2;
1007 edgesWithNodes.elems[0].vtkType = VTK_LINE;
1009 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1010 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1011 edgesWithNodes.elems[1].nbNodes = 2;
1012 edgesWithNodes.elems[1].vtkType = VTK_LINE;
1014 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1015 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1016 edgesWithNodes.elems[2].nbNodes = 2;
1017 edgesWithNodes.elems[2].vtkType = VTK_LINE;
1019 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1020 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1021 edgesWithNodes.elems[3].nbNodes = 2;
1022 edgesWithNodes.elems[3].vtkType = VTK_LINE;
1025 void SMDS_DownQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
1027 //ASSERT((cellId >=0)&& (cellId < _maxId));
1028 //ASSERT(aType == VTK_LINE);
1029 int *faces = &_cellIds[_nbDownCells * cellId];
1030 for (int i = 0; i < _nbDownCells; i++)
1034 faces[i] = lowCellId;
1037 if (faces[i] == lowCellId)
1043 // ---------------------------------------------------------------------------
1045 SMDS_DownQuadQuadrangle::SMDS_DownQuadQuadrangle(SMDS_UnstructuredGrid *grid) :
1046 SMDS_Down2D(grid, 4)
1048 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1049 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1050 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1051 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1055 SMDS_DownQuadQuadrangle::~SMDS_DownQuadQuadrangle()
1059 void SMDS_DownQuadQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
1061 int *nodes = &_tempNodes[_nbNodes * cellId];
1062 edgesWithNodes.nbElems = 4;
1064 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1065 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1066 edgesWithNodes.elems[0].nodeIds[2] = nodes[4];
1067 edgesWithNodes.elems[0].nbNodes = 3;
1068 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
1070 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1071 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1072 edgesWithNodes.elems[1].nodeIds[2] = nodes[5];
1073 edgesWithNodes.elems[1].nbNodes = 3;
1074 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
1076 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1077 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1078 edgesWithNodes.elems[2].nodeIds[2] = nodes[6];
1079 edgesWithNodes.elems[2].nbNodes = 3;
1080 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
1082 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1083 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1084 edgesWithNodes.elems[3].nodeIds[2] = nodes[7];
1085 edgesWithNodes.elems[3].nbNodes = 3;
1086 edgesWithNodes.elems[3].vtkType = VTK_QUADRATIC_EDGE;
1089 void SMDS_DownQuadQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
1091 //ASSERT((cellId >=0)&& (cellId < _maxId));
1092 //ASSERT(aType == VTK_QUADRATIC_EDGE);
1093 int *faces = &_cellIds[_nbDownCells * cellId];
1094 for (int i = 0; i < _nbDownCells; i++)
1098 faces[i] = lowCellId;
1101 if (faces[i] == lowCellId)
1107 // ---------------------------------------------------------------------------
1109 SMDS_DownTetra::SMDS_DownTetra(SMDS_UnstructuredGrid *grid) :
1110 SMDS_Down3D(grid, 4)
1112 _cellTypes.push_back(VTK_TRIANGLE);
1113 _cellTypes.push_back(VTK_TRIANGLE);
1114 _cellTypes.push_back(VTK_TRIANGLE);
1115 _cellTypes.push_back(VTK_TRIANGLE);
1118 SMDS_DownTetra::~SMDS_DownTetra()
1122 void SMDS_DownTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1126 for ( size_t i = 0; i < orderedNodes.size(); i++ )
1127 setNodes.insert(orderedNodes[i]);
1128 //MESSAGE("cellId = " << cellId);
1131 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1132 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1135 int ids[12] = { 0, 1, 2, 0, 3, 1, 2, 3, 0, 1, 3, 2 };
1136 //int ids[12] = { 2, 1, 0, 1, 3, 0, 0, 3, 2, 2, 3, 1 };
1137 for (int k = 0; k < 4; k++)
1140 for (int i = 0; i < 3; i++)
1141 tofind.insert(nodes[ids[3 * k + i]]);
1142 if (setNodes == tofind)
1144 for (int i = 0; i < 3; i++)
1145 orderedNodes[i] = nodes[ids[3 * k + i]];
1149 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
1150 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1151 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1154 void SMDS_DownTetra::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
1156 //ASSERT((cellId >=0)&& (cellId < _maxId));
1157 //ASSERT(aType == VTK_TRIANGLE);
1158 int *faces = &_cellIds[_nbDownCells * cellId];
1159 for (int i = 0; i < _nbDownCells; i++)
1163 faces[i] = lowCellId;
1166 if (faces[i] == lowCellId)
1172 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1173 * The linear tetrahedron is defined by four points.
1174 * @see vtkTetra.h in Filtering.
1175 * @param cellId volumeId in vtkUnstructuredGrid
1176 * @param facesWithNodes vector of face descriptors to be filled
1178 void SMDS_DownTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1180 // --- find point id's of the volume
1183 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1184 _grid->GetCellPoints(cellId, npts, nodes);
1186 // --- create all the ordered list of node id's for each face
1188 facesWithNodes.nbElems = 4;
1190 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1191 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1192 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1193 facesWithNodes.elems[0].nbNodes = 3;
1194 facesWithNodes.elems[0].vtkType = VTK_TRIANGLE;
1196 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1197 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1198 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1199 facesWithNodes.elems[1].nbNodes = 3;
1200 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1202 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1203 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1204 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1205 facesWithNodes.elems[2].nbNodes = 3;
1206 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1208 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1209 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1210 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1211 facesWithNodes.elems[3].nbNodes = 3;
1212 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1215 // ---------------------------------------------------------------------------
1217 SMDS_DownQuadTetra::SMDS_DownQuadTetra(SMDS_UnstructuredGrid *grid) :
1218 SMDS_Down3D(grid, 4)
1220 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1221 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1222 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1223 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1226 SMDS_DownQuadTetra::~SMDS_DownQuadTetra()
1230 void SMDS_DownQuadTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1234 for ( size_t i = 0; i < orderedNodes.size(); i++ )
1235 setNodes.insert(orderedNodes[i]);
1236 //MESSAGE("cellId = " << cellId);
1239 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1240 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1243 int ids[24] = { 0, 1, 2, 4, 5, 6, 0, 3, 1, 7, 8, 4, 2, 3, 0, 9, 7, 6, 1, 3, 2, 8, 9, 5 };
1244 //int ids[24] = { 2, 1, 0, 5, 4, 6, 1, 3, 0, 8, 7, 4, 0, 3, 2, 7, 9, 6, 2, 3, 1, 9, 8, 5 };
1245 for (int k = 0; k < 4; k++)
1248 for (int i = 0; i < 6; i++)
1249 tofind.insert(nodes[ids[6 * k + i]]);
1250 if (setNodes == tofind)
1252 for (int i = 0; i < 6; i++)
1253 orderedNodes[i] = nodes[ids[6 * k + i]];
1257 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
1258 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1259 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1262 void SMDS_DownQuadTetra::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
1264 //ASSERT((cellId >=0)&& (cellId < _maxId));
1265 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1266 int *faces = &_cellIds[_nbDownCells * cellId];
1267 for (int i = 0; i < _nbDownCells; i++)
1271 faces[i] = lowCellId;
1274 if (faces[i] == lowCellId)
1280 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1281 * The ordering of the ten points defining the quadratic tetrahedron cell is point id's (0-3,4-9)
1282 * where id's 0-3 are the four tetrahedron vertices;
1283 * and point id's 4-9 are the mid-edge nodes between (0,1), (1,2), (2,0), (0,3), (1,3), and (2,3).
1284 * @see vtkQuadraticTetra.h in Filtering.
1285 * @param cellId volumeId in vtkUnstructuredGrid
1286 * @param facesWithNodes vector of face descriptors to be filled
1288 void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1290 // --- find point id's of the volume
1293 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1294 _grid->GetCellPoints(cellId, npts, nodes);
1296 // --- create all the ordered list of node id's for each face
1298 facesWithNodes.nbElems = 4;
1300 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1301 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1302 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1303 facesWithNodes.elems[0].nodeIds[3] = nodes[4];
1304 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1305 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1306 facesWithNodes.elems[0].nbNodes = 6;
1307 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_TRIANGLE;
1309 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1310 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1311 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1312 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1313 facesWithNodes.elems[1].nodeIds[4] = nodes[8];
1314 facesWithNodes.elems[1].nodeIds[5] = nodes[7];
1315 facesWithNodes.elems[1].nbNodes = 6;
1316 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1318 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1319 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1320 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1321 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1322 facesWithNodes.elems[2].nodeIds[4] = nodes[9];
1323 facesWithNodes.elems[2].nodeIds[5] = nodes[7];
1324 facesWithNodes.elems[2].nbNodes = 6;
1325 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1327 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1328 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1329 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1330 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
1331 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
1332 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1333 facesWithNodes.elems[3].nbNodes = 6;
1334 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1337 // ---------------------------------------------------------------------------
1339 SMDS_DownPyramid::SMDS_DownPyramid(SMDS_UnstructuredGrid *grid) :
1340 SMDS_Down3D(grid, 5)
1342 _cellTypes.push_back(VTK_QUAD);
1343 _cellTypes.push_back(VTK_TRIANGLE);
1344 _cellTypes.push_back(VTK_TRIANGLE);
1345 _cellTypes.push_back(VTK_TRIANGLE);
1346 _cellTypes.push_back(VTK_TRIANGLE);
1349 SMDS_DownPyramid::~SMDS_DownPyramid()
1353 void SMDS_DownPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1357 for ( size_t i = 0; i < orderedNodes.size(); i++ )
1358 setNodes.insert(orderedNodes[i]);
1359 //MESSAGE("cellId = " << cellId);
1362 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1363 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1366 int ids[16] = { 0, 1, 2, 3, 0, 3, 4, 3, 2, 4, 2, 1, 4, 1, 0, 4 };
1368 // Quadrangular face
1370 for (int i = 0; i < 4; i++)
1371 tofind.insert(nodes[ids[i]]);
1372 if (setNodes == tofind)
1374 for (int i = 0; i < 4; i++)
1375 orderedNodes[i] = nodes[ids[i]];
1379 for (int k = 0; k < 4; k++)
1382 for (int i = 0; i < 3; i++)
1383 tofind.insert(nodes[ids[4 + 3 * k + i]]);
1384 if (setNodes == tofind)
1386 for (int i = 0; i < 3; i++)
1387 orderedNodes[i] = nodes[ids[4 + 3 * k + i]];
1391 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
1392 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1393 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1396 void SMDS_DownPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1398 //ASSERT((cellId >=0) && (cellId < _maxId));
1399 int *faces = &_cellIds[_nbDownCells * cellId];
1400 if (aType == VTK_QUAD)
1404 faces[0] = lowCellId;
1407 if (faces[0] == lowCellId)
1412 //ASSERT(aType == VTK_TRIANGLE);
1413 for (int i = 1; i < _nbDownCells; i++)
1417 faces[i] = lowCellId;
1420 if (faces[i] == lowCellId)
1427 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1428 * The pyramid is defined by the five points (0-4) where (0,1,2,3) is the base of the pyramid which,
1429 * using the right hand rule, forms a quadrilateral whose normal points in the direction of the
1430 * pyramid apex at vertex #4.
1431 * @see vtkPyramid.h in Filtering.
1432 * @param cellId volumeId in vtkUnstructuredGrid
1433 * @param facesWithNodes vector of face descriptors to be filled
1435 void SMDS_DownPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1437 // --- find point id's of the volume
1440 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1441 _grid->GetCellPoints(cellId, npts, nodes);
1443 // --- create all the ordered list of node id's for each face
1445 facesWithNodes.nbElems = 5;
1447 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1448 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1449 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1450 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1451 facesWithNodes.elems[0].nbNodes = 4;
1452 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1454 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1455 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1456 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1457 facesWithNodes.elems[1].nbNodes = 3;
1458 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1460 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1461 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1462 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1463 facesWithNodes.elems[2].nbNodes = 3;
1464 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1466 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1467 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1468 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1469 facesWithNodes.elems[3].nbNodes = 3;
1470 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1472 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1473 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1474 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1475 facesWithNodes.elems[4].nbNodes = 3;
1476 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1479 // ---------------------------------------------------------------------------
1481 SMDS_DownQuadPyramid::SMDS_DownQuadPyramid(SMDS_UnstructuredGrid *grid) :
1482 SMDS_Down3D(grid, 5)
1484 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1485 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1486 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1487 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1488 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1491 SMDS_DownQuadPyramid::~SMDS_DownQuadPyramid()
1495 void SMDS_DownQuadPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1497 // MESSAGE("SMDS_DownQuadPyramid::getOrderedNodesOfFace cellId = " << cellId);
1500 for ( size_t i = 0; i < orderedNodes.size(); i++)
1501 setNodes.insert(orderedNodes[i]);
1502 //MESSAGE("cellId = " << cellId);
1505 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1506 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1509 int ids[32] = { 0, 1, 2, 3, 5, 6, 7, 8,
1510 0, 3, 4, 8, 12, 9, 3, 2, 4, 7 , 11, 12, 2, 1, 4, 6, 10, 11, 1, 0, 4, 5, 9, 10 };
1512 // Quadrangular face
1514 for (int i = 0; i < 8; i++)
1515 tofind.insert(nodes[ids[i]]);
1516 if (setNodes == tofind)
1518 for (int i = 0; i < 8; i++)
1519 orderedNodes[i] = nodes[ids[i]];
1523 for (int k = 0; k < 4; k++)
1526 for (int i = 0; i < 6; i++)
1527 tofind.insert(nodes[ids[8 + 6 * k + i]]);
1528 if (setNodes == tofind)
1530 for (int i = 0; i < 6; i++)
1531 orderedNodes[i] = nodes[ids[8 + 6 * k + i]];
1535 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
1536 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
1537 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1540 void SMDS_DownQuadPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1542 //ASSERT((cellId >=0) && (cellId < _maxId));
1543 int *faces = &_cellIds[_nbDownCells * cellId];
1544 if (aType == VTK_QUADRATIC_QUAD)
1548 faces[0] = lowCellId;
1551 if (faces[0] == lowCellId)
1556 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1557 for (int i = 1; i < _nbDownCells; i++)
1561 faces[i] = lowCellId;
1564 if (faces[i] == lowCellId)
1571 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1572 * The ordering of the thirteen points defining the quadratic pyramid cell is point id's (0-4,5-12)
1573 * where point id's 0-4 are the five corner vertices of the pyramid; followed
1574 * by eight mid-edge nodes (5-12). Note that these mid-edge nodes lie on the edges defined by
1575 * 5(0,1), 6(1,2), 7(2,3), 8(3,0), 9(0,4), 10(1,4), 11(2,4), 12(3,4).
1576 * @see vtkQuadraticPyramid.h in Filtering.
1577 * @param cellId volumeId in vtkUnstructuredGrid
1578 * @param facesWithNodes vector of face descriptors to be filled
1580 void SMDS_DownQuadPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1582 // --- find point id's of the volume
1585 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1586 _grid->GetCellPoints(cellId, npts, nodes);
1588 // --- create all the ordered list of node id's for each face
1590 facesWithNodes.nbElems = 5;
1592 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1593 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1594 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1595 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1596 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1597 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1598 facesWithNodes.elems[0].nodeIds[6] = nodes[7];
1599 facesWithNodes.elems[0].nodeIds[7] = nodes[8];
1600 facesWithNodes.elems[0].nbNodes = 8;
1601 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1603 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1604 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1605 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1606 facesWithNodes.elems[1].nodeIds[3] = nodes[5];
1607 facesWithNodes.elems[1].nodeIds[4] = nodes[10];
1608 facesWithNodes.elems[1].nodeIds[5] = nodes[9];
1609 facesWithNodes.elems[1].nbNodes = 6;
1610 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1612 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1613 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1614 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1615 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1616 facesWithNodes.elems[2].nodeIds[4] = nodes[11];
1617 facesWithNodes.elems[2].nodeIds[5] = nodes[10];
1618 facesWithNodes.elems[2].nbNodes = 6;
1619 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1621 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1622 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1623 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1624 facesWithNodes.elems[3].nodeIds[3] = nodes[7];
1625 facesWithNodes.elems[3].nodeIds[4] = nodes[12];
1626 facesWithNodes.elems[3].nodeIds[5] = nodes[11];
1627 facesWithNodes.elems[3].nbNodes = 6;
1628 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1630 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1631 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1632 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1633 facesWithNodes.elems[4].nodeIds[3] = nodes[8];
1634 facesWithNodes.elems[4].nodeIds[4] = nodes[9];
1635 facesWithNodes.elems[4].nodeIds[5] = nodes[12];
1636 facesWithNodes.elems[4].nbNodes = 6;
1637 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1640 // ---------------------------------------------------------------------------
1642 SMDS_DownPenta::SMDS_DownPenta(SMDS_UnstructuredGrid *grid) :
1643 SMDS_Down3D(grid, 5)
1645 _cellTypes.push_back(VTK_QUAD);
1646 _cellTypes.push_back(VTK_QUAD);
1647 _cellTypes.push_back(VTK_QUAD);
1648 _cellTypes.push_back(VTK_TRIANGLE);
1649 _cellTypes.push_back(VTK_TRIANGLE);
1652 SMDS_DownPenta::~SMDS_DownPenta()
1656 void SMDS_DownPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1660 for ( size_t i = 0; i < orderedNodes.size(); i++)
1661 setNodes.insert(orderedNodes[i]);
1662 //MESSAGE("cellId = " << cellId);
1665 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1666 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1669 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1670 int ids[18] = { 0, 1, 2, 3, 5, 4, 0, 3, 4, 1, 1, 4, 5, 2, 2, 5, 3, 0 };
1673 for (int k = 0; k < 2; k++)
1676 for (int i = 0; i < 3; i++)
1677 tofind.insert(nodes[ids[3 * k + i]]);
1678 if (setNodes == tofind)
1680 for (int i = 0; i < 3; i++)
1681 orderedNodes[i] = nodes[ids[3 * k + i]];
1685 // Quadrangular faces
1686 for (int k = 0; k < 3; k++)
1689 for (int i = 0; i < 4; i++)
1690 tofind.insert(nodes[ids[6 + 4 * k + i]]);
1691 if (setNodes == tofind)
1693 for (int i = 0; i < 4; i++)
1694 orderedNodes[i] = nodes[ids[6 + 4 * k + i]];
1698 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
1699 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1700 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1703 void SMDS_DownPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1705 //ASSERT((cellId >=0) && (cellId < _maxId));
1706 int *faces = &_cellIds[_nbDownCells * cellId];
1707 if (aType == VTK_QUAD)
1708 for (int i = 0; i < 3; i++)
1712 faces[i] = lowCellId;
1715 if (faces[i] == lowCellId)
1720 //ASSERT(aType == VTK_TRIANGLE);
1721 for (int i = 3; i < _nbDownCells; i++)
1725 faces[i] = lowCellId;
1728 if (faces[i] == lowCellId)
1735 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's.
1736 * A wedge or pentahedron consists of two triangular and three quadrilateral faces
1737 * and is defined by the six points (0-5) where (0,1,2) is the base of the wedge which,
1738 * using the right hand rule, forms a triangle whose normal points outward
1739 * (away from the triangular face (3,4,5)).
1740 * @see vtkWedge.h in Filtering
1741 * @param cellId volumeId in vtkUnstructuredGrid
1742 * @param facesWithNodes vector of face descriptors to be filled
1744 void SMDS_DownPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1746 // --- find point id's of the volume
1749 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1750 _grid->GetCellPoints(cellId, npts, nodes);
1752 // --- create all the ordered list of node id's for each face
1754 facesWithNodes.nbElems = 5;
1756 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1757 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1758 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1759 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1760 facesWithNodes.elems[0].nbNodes = 4;
1761 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1763 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1764 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1765 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1766 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1767 facesWithNodes.elems[1].nbNodes = 4;
1768 facesWithNodes.elems[1].vtkType = VTK_QUAD;
1770 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1771 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1772 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1773 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1774 facesWithNodes.elems[2].nbNodes = 4;
1775 facesWithNodes.elems[2].vtkType = VTK_QUAD;
1777 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1778 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1779 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1780 facesWithNodes.elems[3].nbNodes = 3;
1781 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1783 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1784 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1785 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1786 facesWithNodes.elems[4].nbNodes = 3;
1787 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1790 // ---------------------------------------------------------------------------
1792 SMDS_DownQuadPenta::SMDS_DownQuadPenta(SMDS_UnstructuredGrid *grid) :
1793 SMDS_Down3D(grid, 5)
1795 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1796 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1797 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1798 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1799 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1802 SMDS_DownQuadPenta::~SMDS_DownQuadPenta()
1806 void SMDS_DownQuadPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1810 for ( size_t i = 0; i < orderedNodes.size(); i++ )
1811 setNodes.insert(orderedNodes[i]);
1812 //MESSAGE("cellId = " << cellId);
1815 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1816 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1819 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1820 int ids[36] = { 0, 1, 2, 6, 7, 8, 3, 5, 4, 11, 10, 9,
1821 0, 3, 4, 1, 12, 9, 13, 6, 1, 4, 5, 2, 13, 10, 14, 7, 2, 5, 3, 0, 14, 11, 12, 8 };
1824 for (int k = 0; k < 2; k++)
1827 for (int i = 0; i < 6; i++)
1828 tofind.insert(nodes[ids[6 * k + i]]);
1829 if (setNodes == tofind)
1831 for (int i = 0; i < 6; i++)
1832 orderedNodes[i] = nodes[ids[6 * k + i]];
1836 // Quadrangular faces
1837 for (int k = 0; k < 3; k++)
1840 for (int i = 0; i < 8; i++)
1841 tofind.insert(nodes[ids[12 + 8 * k + i]]);
1842 if (setNodes == tofind)
1844 for (int i = 0; i < 8; i++)
1845 orderedNodes[i] = nodes[ids[12 + 8 * k + i]];
1849 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
1850 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1851 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1854 void SMDS_DownQuadPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1856 //ASSERT((cellId >=0) && (cellId < _maxId));
1857 int *faces = &_cellIds[_nbDownCells * cellId];
1858 if (aType == VTK_QUADRATIC_QUAD)
1859 for (int i = 0; i < 3; i++)
1863 faces[i] = lowCellId;
1866 if (faces[i] == lowCellId)
1871 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1872 for (int i = 3; i < _nbDownCells; i++)
1876 faces[i] = lowCellId;
1879 if (faces[i] == lowCellId)
1886 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1887 * The quadratic wedge (or pentahedron) is defined by fifteen points.
1888 * The ordering of the fifteen points defining the cell is point id's (0-5,6-14)
1889 * where point id's 0-5 are the six corner vertices of the wedge, followed by
1890 * nine mid-edge nodes (6-14). Note that these mid-edge nodes lie on the edges defined by
1891 * (0,1), (1,2), (2,0), (3,4), (4,5), (5,3), (0,3), (1,4), (2,5).
1892 * @see vtkQuadraticWedge.h in Filtering
1893 * @param cellId volumeId in vtkUnstructuredGrid
1894 * @param facesWithNodes vector of face descriptors to be filled
1896 void SMDS_DownQuadPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1898 // --- find point id's of the volume
1901 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1902 _grid->GetCellPoints(cellId, npts, nodes);
1904 // --- create all the ordered list of node id's for each face
1906 facesWithNodes.nbElems = 5;
1908 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1909 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1910 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1911 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1912 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
1913 facesWithNodes.elems[0].nodeIds[5] = nodes[14];
1914 facesWithNodes.elems[0].nodeIds[6] = nodes[11];
1915 facesWithNodes.elems[0].nodeIds[7] = nodes[12];
1916 facesWithNodes.elems[0].nbNodes = 8;
1917 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1919 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1920 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1921 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1922 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1923 facesWithNodes.elems[1].nodeIds[4] = nodes[7];
1924 facesWithNodes.elems[1].nodeIds[5] = nodes[14];
1925 facesWithNodes.elems[1].nodeIds[6] = nodes[10];
1926 facesWithNodes.elems[1].nodeIds[7] = nodes[13];
1927 facesWithNodes.elems[1].nbNodes = 8;
1928 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
1930 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1931 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1932 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1933 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1934 facesWithNodes.elems[2].nodeIds[4] = nodes[6];
1935 facesWithNodes.elems[2].nodeIds[5] = nodes[13];
1936 facesWithNodes.elems[2].nodeIds[6] = nodes[9];
1937 facesWithNodes.elems[2].nodeIds[7] = nodes[12];
1938 facesWithNodes.elems[2].nbNodes = 8;
1939 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
1941 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1942 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1943 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1944 facesWithNodes.elems[3].nodeIds[3] = nodes[6];
1945 facesWithNodes.elems[3].nodeIds[4] = nodes[7];
1946 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1947 facesWithNodes.elems[3].nbNodes = 6;
1948 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1950 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1951 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1952 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1953 facesWithNodes.elems[4].nodeIds[3] = nodes[9];
1954 facesWithNodes.elems[4].nodeIds[4] = nodes[10];
1955 facesWithNodes.elems[4].nodeIds[5] = nodes[11];
1956 facesWithNodes.elems[4].nbNodes = 6;
1957 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1960 // ---------------------------------------------------------------------------
1962 SMDS_DownHexa::SMDS_DownHexa(SMDS_UnstructuredGrid *grid) :
1963 SMDS_Down3D(grid, 6)
1965 _cellTypes.push_back(VTK_QUAD);
1966 _cellTypes.push_back(VTK_QUAD);
1967 _cellTypes.push_back(VTK_QUAD);
1968 _cellTypes.push_back(VTK_QUAD);
1969 _cellTypes.push_back(VTK_QUAD);
1970 _cellTypes.push_back(VTK_QUAD);
1973 SMDS_DownHexa::~SMDS_DownHexa()
1977 void SMDS_DownHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1981 for ( size_t i = 0; i < orderedNodes.size(); i++ )
1982 setNodes.insert(orderedNodes[i]);
1983 //MESSAGE("cellId = " << cellId);
1986 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
1987 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1990 //int ids[24] = { 0, 1, 2, 3, 7, 6, 5, 4, 4, 0, 3, 7, 5, 1, 0, 4, 6, 2, 1, 5, 7, 3, 2, 6};
1991 int ids[24] = { 3, 2, 1, 0, 4, 5, 6, 7, 7, 3, 0, 4, 4, 0, 1, 5, 5, 1, 2, 6, 6, 2, 3, 7};
1992 for (int k = 0; k < 6; k++) // loop on the 6 faces
1995 for (int i = 0; i < 4; i++)
1996 tofind.insert(nodes[ids[4 * k + i]]); // node ids of the face i
1997 if (setNodes == tofind)
1999 for (int i = 0; i < 4; i++)
2000 orderedNodes[i] = nodes[ids[4 * k + i]];
2004 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
2005 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2006 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2007 MESSAGE(nodes[4] << " " << nodes[5] << " " << nodes[6] << " " << nodes[7]);
2010 void SMDS_DownHexa::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
2012 //ASSERT((cellId >=0)&& (cellId < _maxId));
2013 int *faces = &_cellIds[_nbDownCells * cellId];
2014 for (int i = 0; i < _nbDownCells; i++)
2018 faces[i] = lowCellId;
2021 if (faces[i] == lowCellId)
2025 // MESSAGE("-------------------------------------> trop de faces ! " << cellId << " " << lowCellId);
2028 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2029 * The hexahedron is defined by the eight points (0-7), where (0,1,2,3) is the base
2030 * of the hexahedron which, using the right hand rule, forms a quadrilateral whose normal
2031 * points in the direction of the opposite face (4,5,6,7).
2032 * @see vtkHexahedron.h in Filtering
2033 * @param cellId volumeId in vtkUnstructuredGrid
2034 * @param facesWithNodes vector of face descriptors to be filled
2036 void SMDS_DownHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2038 // --- find point id's of the volume
2041 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
2042 _grid->GetCellPoints(cellId, npts, nodes);
2044 // --- create all the ordered list of node id's for each face
2046 facesWithNodes.nbElems = 6;
2048 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2049 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2050 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2051 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2052 facesWithNodes.elems[0].nbNodes = 4;
2053 facesWithNodes.elems[0].vtkType = VTK_QUAD;
2055 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2056 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2057 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2058 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2059 facesWithNodes.elems[1].nbNodes = 4;
2060 facesWithNodes.elems[1].vtkType = VTK_QUAD;
2062 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2063 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2064 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2065 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2066 facesWithNodes.elems[2].nbNodes = 4;
2067 facesWithNodes.elems[2].vtkType = VTK_QUAD;
2069 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2070 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2071 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2072 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2073 facesWithNodes.elems[3].nbNodes = 4;
2074 facesWithNodes.elems[3].vtkType = VTK_QUAD;
2076 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2077 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2078 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2079 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2080 facesWithNodes.elems[4].nbNodes = 4;
2081 facesWithNodes.elems[4].vtkType = VTK_QUAD;
2083 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2084 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2085 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2086 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2087 facesWithNodes.elems[5].nbNodes = 4;
2088 facesWithNodes.elems[5].vtkType = VTK_QUAD;
2091 // ---------------------------------------------------------------------------
2093 SMDS_DownQuadHexa::SMDS_DownQuadHexa(SMDS_UnstructuredGrid *grid) :
2094 SMDS_Down3D(grid, 6)
2096 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2097 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2098 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2099 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2100 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2101 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2104 SMDS_DownQuadHexa::~SMDS_DownQuadHexa()
2108 void SMDS_DownQuadHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
2112 for ( size_t i = 0; i < orderedNodes.size(); i++ )
2113 setNodes.insert(orderedNodes[i]);
2114 //MESSAGE("cellId = " << cellId);
2117 vtkIdType const *nodes(nullptr); // will refer to the point id's of the volume
2118 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
2121 //int ids[24] = { 3, 2, 1, 0, 4, 5, 6, 7, 7, 3, 0, 4, 4, 0, 1, 5, 5, 1, 2, 6, 6, 2, 3, 7};
2122 int ids[48] = { 3, 2, 1, 0,10, 9, 8,11, 4, 5, 6, 7,12,13,14,15, 7, 3, 0, 4,19,11,16,15,
2123 4, 0, 1, 5,16, 8,17,12, 5, 1, 2, 6,17, 9,18,13, 6, 2, 3, 7,18,10,19,14};
2124 for (int k = 0; k < 6; k++)
2127 for (int i = 0; i < 8; i++)
2128 tofind.insert(nodes[ids[8 * k + i]]);
2129 if (setNodes == tofind)
2131 for (int i = 0; i < 8; i++)
2132 orderedNodes[i] = nodes[ids[8 * k + i]];
2136 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->FromVtkToSmds(_vtkCellIds[cellId]));
2137 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2138 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2141 void SMDS_DownQuadHexa::addDownCell(int cellId, int lowCellId, unsigned char /*aType*/)
2143 //ASSERT((cellId >=0)&& (cellId < _maxId));
2144 int *faces = &_cellIds[_nbDownCells * cellId];
2145 for (int i = 0; i < _nbDownCells; i++)
2149 faces[i] = lowCellId;
2152 if (faces[i] == lowCellId)
2158 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2159 * The ordering of the twenty points defining the quadratic hexahedron cell is point id's (0-7,8-19)
2160 * where point id's 0-7 are the eight corner vertices of the cube, followed by twelve mid-edge nodes (8-19).
2161 * Note that these mid-edge nodes lie on the edges defined by
2162 * (0,1), (1,2), (2,3), (3,0), (4,5), (5,6), (6,7), (7,4), (0,4), (1,5), (2,6), (3,7).
2163 * @see vtkQuadraticHexahedron.h in Filtering
2164 * @param cellId volumeId in vtkUnstructuredGrid
2165 * @param facesWithNodes vector of face descriptors to be filled
2167 void SMDS_DownQuadHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2169 // --- find point id's of the volume
2172 vtkIdType const *nodes; // will refer to the point id's of the volume
2173 _grid->GetCellPoints(cellId, npts, nodes);
2175 // --- create all the ordered list of node id's for each face
2177 facesWithNodes.nbElems = 6;
2179 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2180 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2181 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2182 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2183 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
2184 facesWithNodes.elems[0].nodeIds[5] = nodes[9];
2185 facesWithNodes.elems[0].nodeIds[6] = nodes[10];
2186 facesWithNodes.elems[0].nodeIds[7] = nodes[11];
2187 facesWithNodes.elems[0].nbNodes = 8;
2188 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
2190 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2191 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2192 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2193 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2194 facesWithNodes.elems[1].nodeIds[4] = nodes[12];
2195 facesWithNodes.elems[1].nodeIds[5] = nodes[13];
2196 facesWithNodes.elems[1].nodeIds[6] = nodes[14];
2197 facesWithNodes.elems[1].nodeIds[7] = nodes[15];
2198 facesWithNodes.elems[1].nbNodes = 8;
2199 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
2201 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2202 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2203 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2204 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2205 facesWithNodes.elems[2].nodeIds[4] = nodes[8];
2206 facesWithNodes.elems[2].nodeIds[5] = nodes[17];
2207 facesWithNodes.elems[2].nodeIds[6] = nodes[12];
2208 facesWithNodes.elems[2].nodeIds[7] = nodes[16];
2209 facesWithNodes.elems[2].nbNodes = 8;
2210 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
2212 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2213 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2214 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2215 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2216 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
2217 facesWithNodes.elems[3].nodeIds[5] = nodes[18];
2218 facesWithNodes.elems[3].nodeIds[6] = nodes[13];
2219 facesWithNodes.elems[3].nodeIds[7] = nodes[17];
2220 facesWithNodes.elems[3].nbNodes = 8;
2221 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_QUAD;
2223 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2224 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2225 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2226 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2227 facesWithNodes.elems[4].nodeIds[4] = nodes[18];
2228 facesWithNodes.elems[4].nodeIds[5] = nodes[14];
2229 facesWithNodes.elems[4].nodeIds[6] = nodes[19];
2230 facesWithNodes.elems[4].nodeIds[7] = nodes[10];
2231 facesWithNodes.elems[4].nbNodes = 8;
2232 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_QUAD;
2234 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2235 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2236 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2237 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2238 facesWithNodes.elems[5].nodeIds[4] = nodes[19];
2239 facesWithNodes.elems[5].nodeIds[5] = nodes[15];
2240 facesWithNodes.elems[5].nodeIds[6] = nodes[16];
2241 facesWithNodes.elems[5].nodeIds[7] = nodes[11];
2242 facesWithNodes.elems[5].nbNodes = 8;
2243 facesWithNodes.elems[5].vtkType = VTK_QUADRATIC_QUAD;
2246 // ---------------------------------------------------------------------------