4 * Created on: Jun 3, 2010
8 #include "SMDS_Downward.hxx"
9 #include "SMDS_Mesh.hxx"
10 #include "utilities.h"
12 #include <vtkCellType.h>
13 #include <vtkCellLinks.h>
19 // ---------------------------------------------------------------------------
21 vector<int> SMDS_Downward::_cellDimension;
23 /*! get the dimension of a cell (1,2,3 for 1D, 2D 3D) given the vtk cell type
25 * @param cellType vtk cell type @see vtkCellType.h
28 int SMDS_Downward::getCellDimension(unsigned char cellType)
30 if (_cellDimension.empty())
32 _cellDimension.resize(VTK_MAXTYPE + 1, 0);
33 _cellDimension[VTK_LINE] = 1;
34 _cellDimension[VTK_QUADRATIC_EDGE] = 1;
35 _cellDimension[VTK_TRIANGLE] = 2;
36 _cellDimension[VTK_QUADRATIC_TRIANGLE] = 2;
37 _cellDimension[VTK_QUAD] = 2;
38 _cellDimension[VTK_QUADRATIC_QUAD] = 2;
39 _cellDimension[VTK_TETRA] = 3;
40 _cellDimension[VTK_QUADRATIC_TETRA] = 3;
41 _cellDimension[VTK_HEXAHEDRON] = 3;
42 _cellDimension[VTK_QUADRATIC_HEXAHEDRON] = 3;
43 _cellDimension[VTK_WEDGE] = 3;
44 _cellDimension[VTK_QUADRATIC_WEDGE] = 3;
45 _cellDimension[VTK_PYRAMID] = 3;
46 _cellDimension[VTK_QUADRATIC_PYRAMID] = 3;
48 return _cellDimension[cellType];
51 // ---------------------------------------------------------------------------
53 /*! Generic constructor for all the downward connectivity structures (one per vtk cell type).
54 * The static structure for cell dimension is set only once.
55 * @param grid unstructured grid associated to the mesh.
56 * @param nbDownCells number of downward entities associated to this vtk type of cell.
59 SMDS_Downward::SMDS_Downward(SMDS_UnstructuredGrid *grid, int nbDownCells) :
60 _grid(grid), _nbDownCells(nbDownCells)
63 this->_cellIds.clear();
64 this->_cellTypes.clear();
65 if (_cellDimension.empty())
67 _cellDimension.resize(VTK_MAXTYPE + 1, 0);
68 _cellDimension[VTK_LINE] = 1;
69 _cellDimension[VTK_QUADRATIC_EDGE] = 1;
70 _cellDimension[VTK_TRIANGLE] = 2;
71 _cellDimension[VTK_QUADRATIC_TRIANGLE] = 2;
72 _cellDimension[VTK_QUAD] = 2;
73 _cellDimension[VTK_QUADRATIC_QUAD] = 2;
74 _cellDimension[VTK_TETRA] = 3;
75 _cellDimension[VTK_QUADRATIC_TETRA] = 3;
76 _cellDimension[VTK_HEXAHEDRON] = 3;
77 _cellDimension[VTK_QUADRATIC_HEXAHEDRON] = 3;
78 _cellDimension[VTK_WEDGE] = 3;
79 _cellDimension[VTK_QUADRATIC_WEDGE] = 3;
80 _cellDimension[VTK_PYRAMID] = 3;
81 _cellDimension[VTK_QUADRATIC_PYRAMID] = 3;
85 SMDS_Downward::~SMDS_Downward()
89 /*! Give or create an entry for downward connectivity structure relative to a cell.
90 * If the entry already exists, just return its id, otherwise, create it.
91 * The internal storage memory is allocated if needed.
92 * The SMDS_UnstructuredGrid::_cellIdToDownId vector is completed for vtkUnstructuredGrid cells.
93 * @param vtkId for a vtkUnstructuredGrid cell or -1 (default) for a created downward cell.
94 * @return the rank in downward[vtkType] structure.
96 int SMDS_Downward::addCell(int vtkId)
100 localId = _grid->CellIdToDownId(vtkId);
104 localId = this->_maxId;
106 this->allocate(_maxId);
109 this->_vtkCellIds[localId] = vtkId;
110 _grid->setCellIdToDownId(vtkId, localId);
112 this->initCell(localId);
116 /*! generic method do nothing. see derived methods
120 void SMDS_Downward::initCell(int cellId)
124 /*! Get the number of downward entities associated to a cell (always the same for a given vtk type of cell)
126 * @param cellId not used here.
129 int SMDS_Downward::getNumberOfDownCells(int cellId)
134 /*! get a pointer on the downward entities id's associated to a cell.
135 * @see SMDS_Downward::getNumberOfDownCells for the number of downward entities.
136 * @see SMDS_Downward::getDownTypes for the vtk cell types associated to the downward entities.
137 * @param cellId index of the cell in the downward structure relative to a given vtk cell type.
138 * @return table of downward entities id's.
140 const int* SMDS_Downward::getDownCells(int cellId)
142 //ASSERT((cellId >=0) && (cellId < _maxId));
143 return &_cellIds[_nbDownCells * cellId];
146 /*! get a list of vtk cell types associated to downward entities of a given cell, in the same order
147 * than the downward entities id's list (@see SMDS_Downward::getDownCells).
149 * @param cellId index of the cell in the downward structure relative to a vtk cell type.
150 * @return table of downward entities types.
152 const unsigned char* SMDS_Downward::getDownTypes(int cellId)
154 return &_cellTypes[0];
157 /*! add a downward entity of dimension n-1 (cell or node) to a given cell.
158 * Actual implementation is done in derived methods.
159 * @param cellId index of the parent cell (dimension n) in the downward structure relative to a vtk cell type.
160 * @param lowCellId index of the children cell to add (dimension n-1)
161 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
163 void SMDS_Downward::addDownCell(int cellId, int lowCellId, unsigned char aType)
165 ASSERT(0); // must be re-implemented in derived class
168 /*! add a downward entity of dimension n+1 to a given cell.
169 * Actual implementation is done in derived methods.
170 * @param cellId index of the children cell (dimension n) in the downward structure relative to a vtk cell type.
171 * @param upCellId index of the parent cell to add (dimension n+1)
172 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
174 void SMDS_Downward::addUpCell(int cellId, int upCellId, unsigned char aType)
176 ASSERT(0); // must be re-implemented in derived class
179 int SMDS_Downward::getNodeSet(int cellId, int* nodeSet)
184 // ---------------------------------------------------------------------------
186 SMDS_Down1D::SMDS_Down1D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
187 SMDS_Downward(grid, nbDownCells)
189 _upCellIdsVector.clear();
190 _upCellTypesVector.clear();
192 _upCellTypes.clear();
193 _upCellIndex.clear();
196 SMDS_Down1D::~SMDS_Down1D()
200 /*! clear vectors used to reference 2D cells containing the edge
204 void SMDS_Down1D::initCell(int cellId)
206 _upCellIdsVector[cellId].clear();
207 _upCellTypesVector[cellId].clear();
210 /*! Resize the downward connectivity storage vector if needed.
212 * @param nbElems total number of elements of the same type required
214 void SMDS_Down1D::allocate(int nbElems)
216 if (nbElems >= _vtkCellIds.size())
218 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
219 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
220 _upCellIdsVector.resize(nbElems + SMDS_Mesh::chunkSize);
221 _upCellTypesVector.resize(nbElems + SMDS_Mesh::chunkSize);
225 void SMDS_Down1D::compactStorage()
227 _cellIds.resize(_nbDownCells * _maxId);
228 _vtkCellIds.resize(_maxId);
231 for (int i = 0; i < _maxId; i++)
232 sizeUpCells += _upCellIdsVector[i].size();
233 _upCellIds.resize(sizeUpCells, -1);
234 _upCellTypes.resize(sizeUpCells);
235 _upCellIndex.resize(_maxId + 1, -1); // id and types of rank i correspond to [ _upCellIndex[i], _upCellIndex[i+1] [
237 for (int i = 0; i < _maxId; i++)
239 _upCellIndex[i] = current;
240 for (int j = 0; j < _upCellIdsVector[i].size(); j++)
242 _upCellIds[current] = _upCellIdsVector[i][j];
243 _upCellTypes[current] = _upCellTypesVector[i][j];
247 _upCellIndex[_maxId] = current;
249 _upCellIdsVector.clear();
250 _upCellTypesVector.clear();
253 void SMDS_Down1D::addUpCell(int cellId, int upCellId, unsigned char aType)
255 //ASSERT((cellId >=0) && (cellId < _maxId));
256 int nbFaces = _upCellIdsVector[cellId].size();
257 for (int i = 0; i < nbFaces; i++)
259 if ((_upCellIdsVector[cellId][i] == upCellId) && (_upCellTypesVector[cellId][i] == aType))
261 return; // already done
264 _upCellIdsVector[cellId].push_back(upCellId);
265 _upCellTypesVector[cellId].push_back(aType);
268 int SMDS_Down1D::getNumberOfUpCells(int cellId)
270 //ASSERT((cellId >=0) && (cellId < _maxId));
271 return _upCellIndex[cellId + 1] - _upCellIndex[cellId];
274 const int* SMDS_Down1D::getUpCells(int cellId)
276 //ASSERT((cellId >=0) && (cellId < _maxId));
277 return &_upCellIds[_upCellIndex[cellId]];
280 const unsigned char* SMDS_Down1D::getUpTypes(int cellId)
282 //ASSERT((cellId >=0) && (cellId < _maxId));
283 return &_upCellTypes[_upCellIndex[cellId]];
286 void SMDS_Down1D::getNodeIds(int cellId, std::set<int>& nodeSet)
288 for (int i = 0; i < _nbDownCells; i++)
289 nodeSet.insert(_cellIds[_nbDownCells * cellId + i]);
292 int SMDS_Down1D::getNodeSet(int cellId, int* nodeSet)
294 for (int i = 0; i < _nbDownCells; i++)
295 nodeSet[i] = _cellIds[_nbDownCells * cellId + i];
299 void SMDS_Down1D::setNodes(int cellId, int vtkId)
302 vtkIdType *pts; // will refer to the point id's of the face
303 _grid->GetCellPoints(vtkId, npts, pts);
304 // MESSAGE(vtkId << " " << npts << " " << _nbDownCells);
305 //ASSERT(npts == _nbDownCells);
306 for (int i = 0; i < npts; i++)
308 _cellIds[_nbDownCells * cellId + i] = pts[i];
312 void SMDS_Down1D::setNodes(int cellId, const int* nodeIds)
314 //ASSERT(nodeIds.size() == _nbDownCells);
315 for (int i = 0; i < _nbDownCells; i++)
317 _cellIds[_nbDownCells * cellId + i] = nodeIds[i];
321 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
322 * We keep in the list the cells that contains all the nodes, we keep only volumes and faces.
323 * @param cellId id of the edge in the downward structure
324 * @param vtkIds vector of vtk id's
325 * @return number of vtk cells (size of vector)
327 int SMDS_Down1D::computeVtkCells(int cellId, std::vector<int>& vtkIds)
331 // --- find all the cells the points belong to, and how many of the points belong to a given cell
333 int *pts = &_cellIds[_nbDownCells * cellId];
334 int ncells = this->computeVtkCells(pts, vtkIds);
338 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
340 * @param pts list of points id's defining an edge
341 * @param vtkIds vector of vtk id's
342 * @return number of vtk cells (size of vector)
344 int SMDS_Down1D::computeVtkCells(int *pts, std::vector<int>& vtkIds)
347 // --- find all the cells the points belong to, and how many of the points belong to a given cell
352 for (int i = 0; i < _nbDownCells; i++)
354 vtkIdType point = pts[i];
355 int numCells = _grid->GetLinks()->GetNcells(point);
356 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
357 for (int j = 0; j < numCells; j++)
359 int vtkCellId = cells[j];
361 for (int k = 0; k < cnt; k++)
363 if (cellIds[k] == vtkCellId)
372 cellIds[cnt] = vtkCellId;
374 // TODO ASSERT(cnt<1000);
380 // --- find the face and volume cells: they contains all the points and are of type volume or face
383 for (int i = 0; i < cnt; i++)
385 if (cellCnt[i] == _nbDownCells)
387 int vtkElemId = cellIds[i];
388 int vtkType = _grid->GetCellType(vtkElemId);
389 if (SMDS_Downward::getCellDimension(vtkType) > 1)
391 vtkIds.push_back(vtkElemId);
400 /*! Build the list of downward faces from a list of vtk cells.
402 * @param cellId id of the edge in the downward structure
403 * @param vtkIds vector of vtk id's
404 * @param downFaces vector of face id's in downward structures
405 * @param downTypes vector of face types
406 * @return number of downward faces
408 int SMDS_Down1D::computeFaces(int cellId, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
410 int *pts = &_cellIds[_nbDownCells * cellId];
411 int nbFaces = this->computeFaces(pts, vtkIds, nbcells, downFaces, downTypes);
415 /*! Build the list of downward faces from a list of vtk cells.
417 * @param pts list of points id's defining an edge
418 * @param vtkIds vector of vtk id's
419 * @param downFaces vector of face id's in downward structures
420 * @param downTypes vector of face types
421 * @return number of downward faces
423 int SMDS_Down1D::computeFaces(int* pts, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
426 for (int i = 0; i < nbcells; i++)
428 int vtkId = vtkIds[i];
429 int vtkType = _grid->GetCellType(vtkId);
430 if (SMDS_Downward::getCellDimension(vtkType) == 2)
432 int faceId = _grid->CellIdToDownId(vtkId);
433 downFaces[cnt] = faceId;
434 downTypes[cnt] = vtkType;
437 else if (SMDS_Downward::getCellDimension(vtkType) == 3)
439 int volId = _grid->CellIdToDownId(vtkId);
440 SMDS_Downward * downvol = _grid->getDownArray(vtkType);
441 //const int *downIds = downvol->getDownCells(volId);
442 const unsigned char* downTypesVol = downvol->getDownTypes(volId);
443 int nbFaces = downvol->getNumberOfDownCells(volId);
444 const int* faceIds = downvol->getDownCells(volId);
445 for (int n = 0; n < nbFaces; n++)
447 SMDS_Down2D *downFace = static_cast<SMDS_Down2D*> (_grid->getDownArray(downTypesVol[n]));
448 bool isInFace = downFace->isInFace(faceIds[n], pts, _nbDownCells);
451 bool alreadySet = false;
452 for (int k = 0; k < cnt; k++)
453 if (faceIds[n] == downFaces[k])
460 downFaces[cnt] = faceIds[n];
461 downTypes[cnt] = downTypesVol[n];
471 // ---------------------------------------------------------------------------
473 SMDS_Down2D::SMDS_Down2D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
474 SMDS_Downward(grid, nbDownCells)
477 _upCellTypes.clear();
482 SMDS_Down2D::~SMDS_Down2D()
486 int SMDS_Down2D::getNumberOfUpCells(int cellId)
489 if (_upCellIds[2 * cellId] >= 0)
491 if (_upCellIds[2 * cellId + 1] >= 0)
496 const int* SMDS_Down2D::getUpCells(int cellId)
498 //ASSERT((cellId >=0) && (cellId < _maxId));
499 return &_upCellIds[2 * cellId];
502 const unsigned char* SMDS_Down2D::getUpTypes(int cellId)
504 //ASSERT((cellId >=0) && (cellId < _maxId));
505 return &_upCellTypes[2 * cellId];
508 void SMDS_Down2D::getNodeIds(int cellId, std::set<int>& nodeSet)
510 for (int i = 0; i < _nbDownCells; i++)
512 int downCellId = _cellIds[_nbDownCells * cellId + i];
513 unsigned char cellType = _cellTypes[i];
514 this->_grid->getDownArray(cellType)->getNodeIds(downCellId, nodeSet);
518 /*! Find in vtkUnstructuredGrid the volumes containing a face already stored in vtkUnstructuredGrid.
519 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
520 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
521 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
522 * @param cellId the face cell id in vkUnstructuredGrid
523 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
524 * @return number of volumes (0, 1 or 2)
526 int SMDS_Down2D::computeVolumeIds(int cellId, int* ids)
528 // --- find point id's of the face
531 vtkIdType *pts; // will refer to the point id's of the face
532 _grid->GetCellPoints(cellId, npts, pts);
534 for (int i = 0; i < npts; i++)
535 nodes.push_back(pts[i]);
536 int nvol = this->computeVolumeIdsFromNodesFace(&nodes[0], npts, ids);
540 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
541 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
542 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
543 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
545 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
546 * @return number of volumes (0, 1 or 2)
548 int SMDS_Down2D::computeVolumeIds(ElemByNodesType& faceByNodes, int* ids)
550 int nvol = this->computeVolumeIdsFromNodesFace(&faceByNodes.nodeIds[0], faceByNodes.nbNodes, ids);
554 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
555 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
556 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
557 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
558 * @param pts array of vtk node id's
559 * @param npts number of nodes
561 * @return number of volumes (0, 1 or 2)
563 int SMDS_Down2D::computeVolumeIdsFromNodesFace(int* pts, int npts, int* ids)
566 // --- find all the cells the points belong to, and how many of the points belong to a given cell
571 for (int i = 0; i < npts; i++)
573 vtkIdType point = pts[i];
574 int numCells = _grid->GetLinks()->GetNcells(point);
575 //MESSAGE("cells pour " << i << " " << numCells);
576 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
577 for (int j = 0; j < numCells; j++)
579 int vtkCellId = cells[j];
581 for (int k = 0; k < cnt; k++)
583 if (cellIds[k] == vtkCellId)
592 cellIds[cnt] = vtkCellId;
594 // TODO ASSERT(cnt<1000);
600 // --- find the volume cells: they contains all the points and are of type volume
603 for (int i = 0; i < cnt; i++)
605 //MESSAGE("cell " << cellIds[i] << " points " << cellCnt[i]);
606 if (cellCnt[i] == npts)
608 int vtkElemId = cellIds[i];
609 int vtkType = _grid->GetCellType(vtkElemId);
610 if (SMDS_Downward::getCellDimension(vtkType) == 3)
612 ids[nvol] = vtkElemId; // store the volume id in given vector
623 void SMDS_Down2D::setTempNodes(int cellId, int vtkId)
626 vtkIdType *pts; // will refer to the point id's of the face
627 _grid->GetCellPoints(vtkId, npts, pts);
628 // MESSAGE(vtkId << " " << npts << " " << _nbNodes);
629 //ASSERT(npts == _nbNodes);
630 for (int i = 0; i < npts; i++)
632 _tempNodes[_nbNodes * cellId + i] = pts[i];
636 void SMDS_Down2D::setTempNodes(int cellId, ElemByNodesType& faceByNodes)
638 for (int i = 0; i < faceByNodes.nbNodes; i++)
639 _tempNodes[_nbNodes * cellId + i] = faceByNodes.nodeIds[i];
642 /*! Find if all the nodes belongs to the face.
644 * @param cellId the face cell Id
645 * @param nodeSet set of node id's to be found in the face list of nodes
648 bool SMDS_Down2D::isInFace(int cellId, int *pts, int npts)
651 int *nodes = &_tempNodes[_nbNodes * cellId];
652 for (int j = 0; j < npts; j++)
655 for (int i = 0; i < _nbNodes; i++)
657 if (nodes[i] == point)
664 return (nbFound == npts);
667 /*! Resize the downward connectivity storage vector if needed.
669 * @param nbElems total number of elements of the same type required
671 void SMDS_Down2D::allocate(int nbElems)
673 if (nbElems >= _vtkCellIds.size())
675 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
676 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
677 _upCellIds.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
678 _upCellTypes.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
679 _tempNodes.resize(_nbNodes * (nbElems + SMDS_Mesh::chunkSize), -1);
683 void SMDS_Down2D::compactStorage()
685 _cellIds.resize(_nbDownCells * _maxId);
686 _upCellIds.resize(2 * _maxId);
687 _upCellTypes.resize(2 * _maxId);
688 _vtkCellIds.resize(_maxId);
692 void SMDS_Down2D::addUpCell(int cellId, int upCellId, unsigned char aType)
694 //ASSERT((cellId >=0)&& (cellId < _maxId));
695 int *vols = &_upCellIds[2 * cellId];
696 unsigned char *types = &_upCellTypes[2 * cellId];
697 for (int i = 0; i < 2; i++)
701 vols[i] = upCellId; // use non affected volume
705 if ((vols[i] == upCellId) && (types[i] == aType)) // already done
711 int SMDS_Down2D::getNodeSet(int cellId, int* nodeSet)
713 for (int i = 0; i < _nbNodes; i++)
714 nodeSet[i] = _tempNodes[_nbNodes * cellId + i];
718 int SMDS_Down2D::FindEdgeByNodes(int cellId, ElemByNodesType& edgeByNodes)
720 int *edges = &_cellIds[_nbDownCells * cellId];
721 for (int i = 0; i < _nbDownCells; i++)
723 if ((edges[i] >= 0) && (edgeByNodes.vtkType == _cellTypes[i]))
726 int npts = this->_grid->getDownArray(edgeByNodes.vtkType)->getNodeSet(edges[i], nodeSet);
728 for (int j = 0; j < npts; j++)
730 int point = edgeByNodes.nodeIds[j];
732 for (int k = 0; k < npts; k++)
734 if (nodeSet[k] == point)
750 // ---------------------------------------------------------------------------
752 SMDS_Down3D::SMDS_Down3D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
753 SMDS_Downward(grid, nbDownCells)
757 SMDS_Down3D::~SMDS_Down3D()
761 void SMDS_Down3D::allocate(int nbElems)
763 if (nbElems >= _vtkCellIds.size())
765 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
766 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
770 void SMDS_Down3D::compactStorage()
772 // nothing to do, size was known before
775 int SMDS_Down3D::getNumberOfUpCells(int cellId)
780 const int* SMDS_Down3D::getUpCells(int cellId)
785 const unsigned char* SMDS_Down3D::getUpTypes(int cellId)
790 void SMDS_Down3D::getNodeIds(int cellId, std::set<int>& nodeSet)
792 int vtkId = this->_vtkCellIds[cellId];
794 vtkIdType *nodes; // will refer to the point id's of the volume
795 _grid->GetCellPoints(vtkId, npts, nodes);
796 for (int i = 0; i < npts; i++)
797 nodeSet.insert(nodes[i]);
800 int SMDS_Down3D::FindFaceByNodes(int cellId, ElemByNodesType& faceByNodes)
802 int *faces = &_cellIds[_nbDownCells * cellId];
805 for (int i = 0; i < _nbDownCells; i++)
807 if ((faces[i] >= 0) && (faceByNodes.vtkType == _cellTypes[i]))
810 npoints = faceByNodes.nbNodes;
813 int npts = this->_grid->getDownArray(faceByNodes.vtkType)->getNodeSet(faces[i], nodeSet);
815 continue; // skip this face
817 for (int j = 0; j < npts; j++)
819 int point = faceByNodes.nodeIds[j];
821 for (int k = 0; k < npts; k++)
823 if (nodeSet[k] == point)
826 break; // point j is in the 2 faces, skip remaining k values
830 break; // point j is not in the 2 faces, skip the remaining tests
839 // ---------------------------------------------------------------------------
841 SMDS_DownEdge::SMDS_DownEdge(SMDS_UnstructuredGrid *grid) :
844 _cellTypes.push_back(VTK_VERTEX);
845 _cellTypes.push_back(VTK_VERTEX);
848 SMDS_DownEdge::~SMDS_DownEdge()
852 // ---------------------------------------------------------------------------
854 SMDS_DownQuadEdge::SMDS_DownQuadEdge(SMDS_UnstructuredGrid *grid) :
857 _cellTypes.push_back(VTK_VERTEX);
858 _cellTypes.push_back(VTK_VERTEX);
859 _cellTypes.push_back(VTK_VERTEX);
862 SMDS_DownQuadEdge::~SMDS_DownQuadEdge()
866 // ---------------------------------------------------------------------------
868 SMDS_DownTriangle::SMDS_DownTriangle(SMDS_UnstructuredGrid *grid) :
871 _cellTypes.push_back(VTK_LINE);
872 _cellTypes.push_back(VTK_LINE);
873 _cellTypes.push_back(VTK_LINE);
877 SMDS_DownTriangle::~SMDS_DownTriangle()
881 void SMDS_DownTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
883 int *nodes = &_tempNodes[_nbNodes * cellId];
884 edgesWithNodes.nbElems = 3;
886 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
887 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
888 edgesWithNodes.elems[0].nbNodes = 2;
889 edgesWithNodes.elems[0].vtkType = VTK_LINE;
891 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
892 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
893 edgesWithNodes.elems[1].nbNodes = 2;
894 edgesWithNodes.elems[1].vtkType = VTK_LINE;
896 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
897 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
898 edgesWithNodes.elems[2].nbNodes = 2;
899 edgesWithNodes.elems[2].vtkType = VTK_LINE;
902 void SMDS_DownTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
904 //ASSERT((cellId >=0)&& (cellId < _maxId));
905 //ASSERT(aType == VTK_LINE);
906 int *faces = &_cellIds[_nbDownCells * cellId];
907 for (int i = 0; i < _nbDownCells; i++)
911 faces[i] = lowCellId;
914 if (faces[i] == lowCellId)
920 // ---------------------------------------------------------------------------
922 SMDS_DownQuadTriangle::SMDS_DownQuadTriangle(SMDS_UnstructuredGrid *grid) :
925 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
926 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
927 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
931 SMDS_DownQuadTriangle::~SMDS_DownQuadTriangle()
935 void SMDS_DownQuadTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
937 int *nodes = &_tempNodes[_nbNodes * cellId];
938 edgesWithNodes.nbElems = 3;
940 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
941 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
942 edgesWithNodes.elems[0].nodeIds[2] = nodes[3];
943 edgesWithNodes.elems[0].nbNodes = 3;
944 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
946 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
947 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
948 edgesWithNodes.elems[1].nodeIds[2] = nodes[4];
949 edgesWithNodes.elems[1].nbNodes = 3;
950 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
952 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
953 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
954 edgesWithNodes.elems[2].nodeIds[2] = nodes[5];
955 edgesWithNodes.elems[2].nbNodes = 3;
956 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
959 void SMDS_DownQuadTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
961 //ASSERT((cellId >=0)&& (cellId < _maxId));
962 //ASSERT(aType == VTK_QUADRATIC_EDGE);
963 int *edges = &_cellIds[_nbDownCells * cellId];
964 for (int i = 0; i < _nbDownCells; i++)
968 edges[i] = lowCellId;
971 if (edges[i] == lowCellId)
977 // ---------------------------------------------------------------------------
979 SMDS_DownQuadrangle::SMDS_DownQuadrangle(SMDS_UnstructuredGrid *grid) :
982 _cellTypes.push_back(VTK_LINE);
983 _cellTypes.push_back(VTK_LINE);
984 _cellTypes.push_back(VTK_LINE);
985 _cellTypes.push_back(VTK_LINE);
989 SMDS_DownQuadrangle::~SMDS_DownQuadrangle()
993 void SMDS_DownQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
995 int *nodes = &_tempNodes[_nbNodes * cellId];
996 edgesWithNodes.nbElems = 4;
998 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
999 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1000 edgesWithNodes.elems[0].nbNodes = 2;
1001 edgesWithNodes.elems[0].vtkType = VTK_LINE;
1003 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1004 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1005 edgesWithNodes.elems[1].nbNodes = 2;
1006 edgesWithNodes.elems[1].vtkType = VTK_LINE;
1008 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1009 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1010 edgesWithNodes.elems[2].nbNodes = 2;
1011 edgesWithNodes.elems[2].vtkType = VTK_LINE;
1013 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1014 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1015 edgesWithNodes.elems[3].nbNodes = 2;
1016 edgesWithNodes.elems[3].vtkType = VTK_LINE;
1019 void SMDS_DownQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1021 //ASSERT((cellId >=0)&& (cellId < _maxId));
1022 //ASSERT(aType == VTK_LINE);
1023 int *faces = &_cellIds[_nbDownCells * cellId];
1024 for (int i = 0; i < _nbDownCells; i++)
1028 faces[i] = lowCellId;
1031 if (faces[i] == lowCellId)
1037 // ---------------------------------------------------------------------------
1039 SMDS_DownQuadQuadrangle::SMDS_DownQuadQuadrangle(SMDS_UnstructuredGrid *grid) :
1040 SMDS_Down2D(grid, 4)
1042 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1043 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1044 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1045 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1049 SMDS_DownQuadQuadrangle::~SMDS_DownQuadQuadrangle()
1053 void SMDS_DownQuadQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
1055 int *nodes = &_tempNodes[_nbNodes * cellId];
1056 edgesWithNodes.nbElems = 4;
1058 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1059 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1060 edgesWithNodes.elems[0].nodeIds[2] = nodes[4];
1061 edgesWithNodes.elems[0].nbNodes = 3;
1062 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
1064 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1065 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1066 edgesWithNodes.elems[1].nodeIds[2] = nodes[5];
1067 edgesWithNodes.elems[1].nbNodes = 3;
1068 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
1070 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1071 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1072 edgesWithNodes.elems[2].nodeIds[2] = nodes[6];
1073 edgesWithNodes.elems[2].nbNodes = 3;
1074 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
1076 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1077 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1078 edgesWithNodes.elems[3].nodeIds[2] = nodes[7];
1079 edgesWithNodes.elems[3].nbNodes = 3;
1080 edgesWithNodes.elems[3].vtkType = VTK_QUADRATIC_EDGE;
1083 void SMDS_DownQuadQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1085 //ASSERT((cellId >=0)&& (cellId < _maxId));
1086 //ASSERT(aType == VTK_QUADRATIC_EDGE);
1087 int *faces = &_cellIds[_nbDownCells * cellId];
1088 for (int i = 0; i < _nbDownCells; i++)
1092 faces[i] = lowCellId;
1095 if (faces[i] == lowCellId)
1101 // ---------------------------------------------------------------------------
1103 SMDS_DownTetra::SMDS_DownTetra(SMDS_UnstructuredGrid *grid) :
1104 SMDS_Down3D(grid, 4)
1106 _cellTypes.push_back(VTK_TRIANGLE);
1107 _cellTypes.push_back(VTK_TRIANGLE);
1108 _cellTypes.push_back(VTK_TRIANGLE);
1109 _cellTypes.push_back(VTK_TRIANGLE);
1112 SMDS_DownTetra::~SMDS_DownTetra()
1116 void SMDS_DownTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1120 for (int i = 0; i < orderedNodes.size(); i++)
1121 setNodes.insert(orderedNodes[i]);
1122 //MESSAGE("cellId = " << cellId);
1125 vtkIdType *nodes; // will refer to the point id's of the volume
1126 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1129 int ids[12] = { 0, 1, 2, 0, 3, 1, 2, 3, 0, 1, 3, 2 };
1130 //int ids[12] = { 2, 1, 0, 1, 3, 0, 0, 3, 2, 2, 3, 1 };
1131 for (int k = 0; k < 4; k++)
1134 for (int i = 0; i < 3; i++)
1135 tofind.insert(nodes[ids[3 * k + i]]);
1136 if (setNodes == tofind)
1138 for (int i = 0; i < 3; i++)
1139 orderedNodes[i] = nodes[ids[3 * k + i]];
1143 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1144 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1145 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1148 void SMDS_DownTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1150 //ASSERT((cellId >=0)&& (cellId < _maxId));
1151 //ASSERT(aType == VTK_TRIANGLE);
1152 int *faces = &_cellIds[_nbDownCells * cellId];
1153 for (int i = 0; i < _nbDownCells; i++)
1157 faces[i] = lowCellId;
1160 if (faces[i] == lowCellId)
1166 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1167 * The linear tetrahedron is defined by four points.
1168 * @see vtkTetra.h in Filtering.
1169 * @param cellId volumeId in vtkUnstructuredGrid
1170 * @param facesWithNodes vector of face descriptors to be filled
1172 void SMDS_DownTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1174 // --- find point id's of the volume
1177 vtkIdType *nodes; // will refer to the point id's of the volume
1178 _grid->GetCellPoints(cellId, npts, nodes);
1180 // --- create all the ordered list of node id's for each face
1182 facesWithNodes.nbElems = 4;
1184 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1185 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1186 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1187 facesWithNodes.elems[0].nbNodes = 3;
1188 facesWithNodes.elems[0].vtkType = VTK_TRIANGLE;
1190 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1191 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1192 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1193 facesWithNodes.elems[1].nbNodes = 3;
1194 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1196 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1197 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1198 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1199 facesWithNodes.elems[2].nbNodes = 3;
1200 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1202 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1203 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1204 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1205 facesWithNodes.elems[3].nbNodes = 3;
1206 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1209 // ---------------------------------------------------------------------------
1211 SMDS_DownQuadTetra::SMDS_DownQuadTetra(SMDS_UnstructuredGrid *grid) :
1212 SMDS_Down3D(grid, 4)
1214 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1215 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1216 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1217 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1220 SMDS_DownQuadTetra::~SMDS_DownQuadTetra()
1224 void SMDS_DownQuadTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1228 for (int i = 0; i < orderedNodes.size(); i++)
1229 setNodes.insert(orderedNodes[i]);
1230 //MESSAGE("cellId = " << cellId);
1233 vtkIdType *nodes; // will refer to the point id's of the volume
1234 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1237 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 };
1238 //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 };
1239 for (int k = 0; k < 4; k++)
1242 for (int i = 0; i < 6; i++)
1243 tofind.insert(nodes[ids[6 * k + i]]);
1244 if (setNodes == tofind)
1246 for (int i = 0; i < 6; i++)
1247 orderedNodes[i] = nodes[ids[6 * k + i]];
1251 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1252 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1253 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1256 void SMDS_DownQuadTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1258 //ASSERT((cellId >=0)&& (cellId < _maxId));
1259 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1260 int *faces = &_cellIds[_nbDownCells * cellId];
1261 for (int i = 0; i < _nbDownCells; i++)
1265 faces[i] = lowCellId;
1268 if (faces[i] == lowCellId)
1274 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1275 * The ordering of the ten points defining the quadratic tetrahedron cell is point id's (0-3,4-9)
1276 * where id's 0-3 are the four tetrahedron vertices;
1277 * 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).
1278 * @see vtkQuadraticTetra.h in Filtering.
1279 * @param cellId volumeId in vtkUnstructuredGrid
1280 * @param facesWithNodes vector of face descriptors to be filled
1282 void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1284 // --- find point id's of the volume
1287 vtkIdType *nodes; // will refer to the point id's of the volume
1288 _grid->GetCellPoints(cellId, npts, nodes);
1290 // --- create all the ordered list of node id's for each face
1292 facesWithNodes.nbElems = 4;
1294 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1295 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1296 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1297 facesWithNodes.elems[0].nodeIds[3] = nodes[4];
1298 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1299 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1300 facesWithNodes.elems[0].nbNodes = 6;
1301 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_TRIANGLE;
1303 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1304 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1305 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1306 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1307 facesWithNodes.elems[1].nodeIds[4] = nodes[8];
1308 facesWithNodes.elems[1].nodeIds[5] = nodes[7];
1309 facesWithNodes.elems[1].nbNodes = 6;
1310 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1312 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1313 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1314 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1315 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1316 facesWithNodes.elems[2].nodeIds[4] = nodes[9];
1317 facesWithNodes.elems[2].nodeIds[5] = nodes[7];
1318 facesWithNodes.elems[2].nbNodes = 6;
1319 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1321 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1322 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1323 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1324 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
1325 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
1326 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1327 facesWithNodes.elems[3].nbNodes = 6;
1328 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1331 // ---------------------------------------------------------------------------
1333 SMDS_DownPyramid::SMDS_DownPyramid(SMDS_UnstructuredGrid *grid) :
1334 SMDS_Down3D(grid, 5)
1336 _cellTypes.push_back(VTK_QUAD);
1337 _cellTypes.push_back(VTK_TRIANGLE);
1338 _cellTypes.push_back(VTK_TRIANGLE);
1339 _cellTypes.push_back(VTK_TRIANGLE);
1340 _cellTypes.push_back(VTK_TRIANGLE);
1343 SMDS_DownPyramid::~SMDS_DownPyramid()
1347 void SMDS_DownPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1351 for (int i = 0; i < orderedNodes.size(); i++)
1352 setNodes.insert(orderedNodes[i]);
1353 //MESSAGE("cellId = " << cellId);
1356 vtkIdType *nodes; // will refer to the point id's of the volume
1357 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1360 int ids[16] = { 0, 1, 2, 3, 0, 3, 4, 3, 2, 4, 2, 1, 4, 1, 0, 4 };
1363 for (int i = 0; i < 4; i++)
1364 tofind.insert(nodes[ids[i]]);
1365 if (setNodes == tofind)
1367 for (int i = 0; i < 4; i++)
1368 orderedNodes[i] = nodes[ids[i]];
1371 for (int k = 0; k < 4; k++)
1374 for (int i = 0; i < 3; i++)
1375 tofind.insert(nodes[ids[4 + 3 * k + i]]);
1376 if (setNodes == tofind)
1378 for (int i = 0; i < 3; i++)
1379 orderedNodes[i] = nodes[ids[4 + 3 * k + i]];
1383 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1384 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1385 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1388 void SMDS_DownPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1390 //ASSERT((cellId >=0) && (cellId < _maxId));
1391 int *faces = &_cellIds[_nbDownCells * cellId];
1392 if (aType == VTK_QUAD)
1396 faces[0] = lowCellId;
1399 if (faces[0] == lowCellId)
1404 //ASSERT(aType == VTK_TRIANGLE);
1405 for (int i = 1; i < _nbDownCells; i++)
1409 faces[i] = lowCellId;
1412 if (faces[i] == lowCellId)
1419 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1420 * The pyramid is defined by the five points (0-4) where (0,1,2,3) is the base of the pyramid which,
1421 * using the right hand rule, forms a quadrilateral whose normal points in the direction of the
1422 * pyramid apex at vertex #4.
1423 * @see vtkPyramid.h in Filtering.
1424 * @param cellId volumeId in vtkUnstructuredGrid
1425 * @param facesWithNodes vector of face descriptors to be filled
1427 void SMDS_DownPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1429 // --- find point id's of the volume
1432 vtkIdType *nodes; // will refer to the point id's of the volume
1433 _grid->GetCellPoints(cellId, npts, nodes);
1435 // --- create all the ordered list of node id's for each face
1437 facesWithNodes.nbElems = 5;
1439 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1440 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1441 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1442 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1443 facesWithNodes.elems[0].nbNodes = 4;
1444 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1446 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1447 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1448 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1449 facesWithNodes.elems[1].nbNodes = 3;
1450 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1452 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1453 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1454 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1455 facesWithNodes.elems[2].nbNodes = 3;
1456 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1458 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1459 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1460 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1461 facesWithNodes.elems[3].nbNodes = 3;
1462 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1464 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1465 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1466 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1467 facesWithNodes.elems[4].nbNodes = 3;
1468 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1471 // ---------------------------------------------------------------------------
1473 SMDS_DownQuadPyramid::SMDS_DownQuadPyramid(SMDS_UnstructuredGrid *grid) :
1474 SMDS_Down3D(grid, 5)
1476 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1477 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1478 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1479 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1480 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1483 SMDS_DownQuadPyramid::~SMDS_DownQuadPyramid()
1487 void SMDS_DownQuadPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1491 for (int i = 0; i < orderedNodes.size(); i++)
1492 setNodes.insert(orderedNodes[i]);
1493 //MESSAGE("cellId = " << cellId);
1496 vtkIdType *nodes; // will refer to the point id's of the volume
1497 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1500 int ids[32] = { 0, 1, 2, 3, 5, 6, 7, 8,
1501 0, 3, 4, 8, 12, 9, 3, 2, 4, 7 , 11, 12, 2, 1, 4, 6, 10, 11, 1, 0, 4, 5, 9, 10 };
1504 for (int i = 0; i < 4; i++)
1505 tofind.insert(nodes[ids[i]]);
1506 if (setNodes == tofind)
1508 for (int i = 0; i < 8; i++)
1509 orderedNodes[i] = nodes[ids[i]];
1512 for (int k = 0; k < 4; k++)
1515 for (int i = 0; i < 6; i++)
1516 tofind.insert(nodes[ids[8 + 6 * k + i]]);
1517 if (setNodes == tofind)
1519 for (int i = 0; i < 6; i++)
1520 orderedNodes[i] = nodes[ids[8 + 6 * k + i]];
1524 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1525 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1526 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1529 void SMDS_DownQuadPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1531 //ASSERT((cellId >=0) && (cellId < _maxId));
1532 int *faces = &_cellIds[_nbDownCells * cellId];
1533 if (aType == VTK_QUADRATIC_QUAD)
1537 faces[0] = lowCellId;
1540 if (faces[0] == lowCellId)
1545 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1546 for (int i = 1; i < _nbDownCells; i++)
1550 faces[i] = lowCellId;
1553 if (faces[i] == lowCellId)
1560 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1561 * The ordering of the thirteen points defining the quadratic pyramid cell is point id's (0-4,5-12)
1562 * where point id's 0-4 are the five corner vertices of the pyramid; followed
1563 * by eight mid-edge nodes (5-12). Note that these mid-edge nodes lie on the edges defined by
1564 * 5(0,1), 6(1,2), 7(2,3), 8(3,0), 9(0,4), 10(1,4), 11(2,4), 12(3,4).
1565 * @see vtkQuadraticPyramid.h in Filtering.
1566 * @param cellId volumeId in vtkUnstructuredGrid
1567 * @param facesWithNodes vector of face descriptors to be filled
1569 void SMDS_DownQuadPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1571 // --- find point id's of the volume
1574 vtkIdType *nodes; // will refer to the point id's of the volume
1575 _grid->GetCellPoints(cellId, npts, nodes);
1577 // --- create all the ordered list of node id's for each face
1579 facesWithNodes.nbElems = 5;
1581 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1582 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1583 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1584 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1585 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1586 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1587 facesWithNodes.elems[0].nodeIds[6] = nodes[7];
1588 facesWithNodes.elems[0].nodeIds[7] = nodes[8];
1589 facesWithNodes.elems[0].nbNodes = 8;
1590 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1592 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1593 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1594 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1595 facesWithNodes.elems[1].nodeIds[3] = nodes[5];
1596 facesWithNodes.elems[1].nodeIds[4] = nodes[10];
1597 facesWithNodes.elems[1].nodeIds[5] = nodes[9];
1598 facesWithNodes.elems[1].nbNodes = 6;
1599 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1601 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1602 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1603 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1604 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1605 facesWithNodes.elems[2].nodeIds[4] = nodes[11];
1606 facesWithNodes.elems[2].nodeIds[5] = nodes[10];
1607 facesWithNodes.elems[2].nbNodes = 6;
1608 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1610 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1611 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1612 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1613 facesWithNodes.elems[3].nodeIds[3] = nodes[7];
1614 facesWithNodes.elems[3].nodeIds[4] = nodes[12];
1615 facesWithNodes.elems[3].nodeIds[5] = nodes[11];
1616 facesWithNodes.elems[3].nbNodes = 6;
1617 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1619 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1620 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1621 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1622 facesWithNodes.elems[4].nodeIds[3] = nodes[8];
1623 facesWithNodes.elems[4].nodeIds[4] = nodes[9];
1624 facesWithNodes.elems[4].nodeIds[5] = nodes[12];
1625 facesWithNodes.elems[4].nbNodes = 6;
1626 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1629 // ---------------------------------------------------------------------------
1631 SMDS_DownPenta::SMDS_DownPenta(SMDS_UnstructuredGrid *grid) :
1632 SMDS_Down3D(grid, 5)
1634 _cellTypes.push_back(VTK_QUAD);
1635 _cellTypes.push_back(VTK_QUAD);
1636 _cellTypes.push_back(VTK_QUAD);
1637 _cellTypes.push_back(VTK_TRIANGLE);
1638 _cellTypes.push_back(VTK_TRIANGLE);
1641 SMDS_DownPenta::~SMDS_DownPenta()
1645 void SMDS_DownPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1649 for (int i = 0; i < orderedNodes.size(); i++)
1650 setNodes.insert(orderedNodes[i]);
1651 //MESSAGE("cellId = " << cellId);
1654 vtkIdType *nodes; // will refer to the point id's of the volume
1655 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1658 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1659 int ids[18] = { 0, 1, 2, 3, 5, 4, 0, 3, 4, 1, 1, 4, 5, 2, 2, 5, 3, 0 };
1661 for (int k = 0; k < 2; k++)
1664 for (int i = 0; i < 3; i++)
1665 tofind.insert(nodes[ids[3 * k + i]]);
1666 if (setNodes == tofind)
1668 for (int i = 0; i < 3; i++)
1669 orderedNodes[i] = nodes[ids[3 * k + i]];
1673 for (int k = 0; k < 3; k++)
1676 for (int i = 0; i < 4; i++)
1677 tofind.insert(nodes[ids[6 + 4 * k + i]]);
1678 if (setNodes == tofind)
1680 for (int i = 0; i < 4; i++)
1681 orderedNodes[i] = nodes[ids[6 + 4 * k + i]];
1685 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1686 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1687 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1690 void SMDS_DownPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1692 //ASSERT((cellId >=0) && (cellId < _maxId));
1693 int *faces = &_cellIds[_nbDownCells * cellId];
1694 if (aType == VTK_QUAD)
1695 for (int i = 0; i < 3; i++)
1699 faces[i] = lowCellId;
1702 if (faces[i] == lowCellId)
1707 //ASSERT(aType == VTK_TRIANGLE);
1708 for (int i = 3; i < _nbDownCells; i++)
1712 faces[i] = lowCellId;
1715 if (faces[i] == lowCellId)
1722 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's.
1723 * A wedge or pentahedron consists of two triangular and three quadrilateral faces
1724 * and is defined by the six points (0-5) where (0,1,2) is the base of the wedge which,
1725 * using the right hand rule, forms a triangle whose normal points outward
1726 * (away from the triangular face (3,4,5)).
1727 * @see vtkWedge.h in Filtering
1728 * @param cellId volumeId in vtkUnstructuredGrid
1729 * @param facesWithNodes vector of face descriptors to be filled
1731 void SMDS_DownPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1733 // --- find point id's of the volume
1736 vtkIdType *nodes; // will refer to the point id's of the volume
1737 _grid->GetCellPoints(cellId, npts, nodes);
1739 // --- create all the ordered list of node id's for each face
1741 facesWithNodes.nbElems = 5;
1743 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1744 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1745 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1746 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1747 facesWithNodes.elems[0].nbNodes = 4;
1748 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1750 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1751 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1752 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1753 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1754 facesWithNodes.elems[1].nbNodes = 4;
1755 facesWithNodes.elems[1].vtkType = VTK_QUAD;
1757 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1758 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1759 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1760 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1761 facesWithNodes.elems[2].nbNodes = 4;
1762 facesWithNodes.elems[2].vtkType = VTK_QUAD;
1764 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1765 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1766 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1767 facesWithNodes.elems[3].nbNodes = 3;
1768 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1770 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1771 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1772 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1773 facesWithNodes.elems[4].nbNodes = 3;
1774 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1777 // ---------------------------------------------------------------------------
1779 SMDS_DownQuadPenta::SMDS_DownQuadPenta(SMDS_UnstructuredGrid *grid) :
1780 SMDS_Down3D(grid, 5)
1782 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1783 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1784 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1785 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1786 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1789 SMDS_DownQuadPenta::~SMDS_DownQuadPenta()
1793 void SMDS_DownQuadPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1797 for (int i = 0; i < orderedNodes.size(); i++)
1798 setNodes.insert(orderedNodes[i]);
1799 //MESSAGE("cellId = " << cellId);
1802 vtkIdType *nodes; // will refer to the point id's of the volume
1803 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1806 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1807 int ids[36] = { 0, 1, 2, 6, 7, 8, 3, 5, 4, 11, 10, 9,
1808 0, 3, 4, 1, 12, 9, 13, 6, 1, 4, 5, 2, 13, 10, 14, 7, 2, 5, 3, 0, 14, 11, 12, 8 };
1810 for (int k = 0; k < 2; k++)
1813 for (int i = 0; i < 6; i++)
1814 tofind.insert(nodes[ids[6 * k + i]]);
1815 if (setNodes == tofind)
1817 for (int i = 0; i < 6; i++)
1818 orderedNodes[i] = nodes[ids[6 * k + i]];
1822 for (int k = 0; k < 3; k++)
1825 for (int i = 0; i < 8; i++)
1826 tofind.insert(nodes[ids[12 + 8 * k + i]]);
1827 if (setNodes == tofind)
1829 for (int i = 0; i < 8; i++)
1830 orderedNodes[i] = nodes[ids[12 + 8 * k + i]];
1834 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1835 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1836 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1839 void SMDS_DownQuadPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1841 //ASSERT((cellId >=0) && (cellId < _maxId));
1842 int *faces = &_cellIds[_nbDownCells * cellId];
1843 if (aType == VTK_QUADRATIC_QUAD)
1844 for (int i = 0; i < 3; i++)
1848 faces[i] = lowCellId;
1851 if (faces[i] == lowCellId)
1856 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1857 for (int i = 3; i < _nbDownCells; i++)
1861 faces[i] = lowCellId;
1864 if (faces[i] == lowCellId)
1871 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1872 * The quadratic wedge (or pentahedron) is defined by fifteen points.
1873 * The ordering of the fifteen points defining the cell is point id's (0-5,6-14)
1874 * where point id's 0-5 are the six corner vertices of the wedge, followed by
1875 * nine mid-edge nodes (6-14). Note that these mid-edge nodes lie on the edges defined by
1876 * (0,1), (1,2), (2,0), (3,4), (4,5), (5,3), (0,3), (1,4), (2,5).
1877 * @see vtkQuadraticWedge.h in Filtering
1878 * @param cellId volumeId in vtkUnstructuredGrid
1879 * @param facesWithNodes vector of face descriptors to be filled
1881 void SMDS_DownQuadPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1883 // --- find point id's of the volume
1886 vtkIdType *nodes; // will refer to the point id's of the volume
1887 _grid->GetCellPoints(cellId, npts, nodes);
1889 // --- create all the ordered list of node id's for each face
1891 facesWithNodes.nbElems = 5;
1893 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1894 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1895 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1896 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1897 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
1898 facesWithNodes.elems[0].nodeIds[5] = nodes[14];
1899 facesWithNodes.elems[0].nodeIds[6] = nodes[11];
1900 facesWithNodes.elems[0].nodeIds[7] = nodes[12];
1901 facesWithNodes.elems[0].nbNodes = 8;
1902 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1904 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1905 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1906 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1907 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1908 facesWithNodes.elems[1].nodeIds[4] = nodes[7];
1909 facesWithNodes.elems[1].nodeIds[5] = nodes[14];
1910 facesWithNodes.elems[1].nodeIds[6] = nodes[10];
1911 facesWithNodes.elems[1].nodeIds[7] = nodes[13];
1912 facesWithNodes.elems[1].nbNodes = 8;
1913 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
1915 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1916 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1917 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1918 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1919 facesWithNodes.elems[2].nodeIds[4] = nodes[6];
1920 facesWithNodes.elems[2].nodeIds[5] = nodes[13];
1921 facesWithNodes.elems[2].nodeIds[6] = nodes[9];
1922 facesWithNodes.elems[2].nodeIds[7] = nodes[12];
1923 facesWithNodes.elems[2].nbNodes = 8;
1924 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
1926 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1927 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1928 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1929 facesWithNodes.elems[3].nodeIds[3] = nodes[6];
1930 facesWithNodes.elems[3].nodeIds[4] = nodes[7];
1931 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1932 facesWithNodes.elems[3].nbNodes = 6;
1933 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1935 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1936 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1937 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1938 facesWithNodes.elems[4].nodeIds[3] = nodes[9];
1939 facesWithNodes.elems[4].nodeIds[4] = nodes[10];
1940 facesWithNodes.elems[4].nodeIds[5] = nodes[11];
1941 facesWithNodes.elems[4].nbNodes = 6;
1942 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1945 // ---------------------------------------------------------------------------
1947 SMDS_DownHexa::SMDS_DownHexa(SMDS_UnstructuredGrid *grid) :
1948 SMDS_Down3D(grid, 6)
1950 _cellTypes.push_back(VTK_QUAD);
1951 _cellTypes.push_back(VTK_QUAD);
1952 _cellTypes.push_back(VTK_QUAD);
1953 _cellTypes.push_back(VTK_QUAD);
1954 _cellTypes.push_back(VTK_QUAD);
1955 _cellTypes.push_back(VTK_QUAD);
1958 SMDS_DownHexa::~SMDS_DownHexa()
1962 void SMDS_DownHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1966 for (int i = 0; i < orderedNodes.size(); i++)
1967 setNodes.insert(orderedNodes[i]);
1968 //MESSAGE("cellId = " << cellId);
1971 vtkIdType *nodes; // will refer to the point id's of the volume
1972 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1975 //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};
1976 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};
1977 for (int k = 0; k < 6; k++) // loop on the 6 faces
1980 for (int i = 0; i < 4; i++)
1981 tofind.insert(nodes[ids[4 * k + i]]); // node ids of the face i
1982 if (setNodes == tofind)
1984 for (int i = 0; i < 4; i++)
1985 orderedNodes[i] = nodes[ids[4 * k + i]];
1989 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1990 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
1991 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1992 MESSAGE(nodes[4] << " " << nodes[5] << " " << nodes[6] << " " << nodes[7]);
1995 void SMDS_DownHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
1997 //ASSERT((cellId >=0)&& (cellId < _maxId));
1998 int *faces = &_cellIds[_nbDownCells * cellId];
1999 for (int i = 0; i < _nbDownCells; i++)
2003 faces[i] = lowCellId;
2006 if (faces[i] == lowCellId)
2010 // MESSAGE("-------------------------------------> trop de faces ! " << cellId << " " << lowCellId);
2013 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2014 * The hexahedron is defined by the eight points (0-7), where (0,1,2,3) is the base
2015 * of the hexahedron which, using the right hand rule, forms a quadrilateral whose normal
2016 * points in the direction of the opposite face (4,5,6,7).
2017 * @see vtkHexahedron.h in Filtering
2018 * @param cellId volumeId in vtkUnstructuredGrid
2019 * @param facesWithNodes vector of face descriptors to be filled
2021 void SMDS_DownHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2023 // --- find point id's of the volume
2026 vtkIdType *nodes; // will refer to the point id's of the volume
2027 _grid->GetCellPoints(cellId, npts, nodes);
2029 // --- create all the ordered list of node id's for each face
2031 facesWithNodes.nbElems = 6;
2033 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2034 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2035 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2036 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2037 facesWithNodes.elems[0].nbNodes = 4;
2038 facesWithNodes.elems[0].vtkType = VTK_QUAD;
2040 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2041 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2042 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2043 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2044 facesWithNodes.elems[1].nbNodes = 4;
2045 facesWithNodes.elems[1].vtkType = VTK_QUAD;
2047 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2048 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2049 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2050 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2051 facesWithNodes.elems[2].nbNodes = 4;
2052 facesWithNodes.elems[2].vtkType = VTK_QUAD;
2054 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2055 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2056 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2057 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2058 facesWithNodes.elems[3].nbNodes = 4;
2059 facesWithNodes.elems[3].vtkType = VTK_QUAD;
2061 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2062 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2063 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2064 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2065 facesWithNodes.elems[4].nbNodes = 4;
2066 facesWithNodes.elems[4].vtkType = VTK_QUAD;
2068 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2069 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2070 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2071 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2072 facesWithNodes.elems[5].nbNodes = 4;
2073 facesWithNodes.elems[5].vtkType = VTK_QUAD;
2076 // ---------------------------------------------------------------------------
2078 SMDS_DownQuadHexa::SMDS_DownQuadHexa(SMDS_UnstructuredGrid *grid) :
2079 SMDS_Down3D(grid, 6)
2081 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2082 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2083 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2084 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2085 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2086 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2089 SMDS_DownQuadHexa::~SMDS_DownQuadHexa()
2093 void SMDS_DownQuadHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
2097 for (int i = 0; i < orderedNodes.size(); i++)
2098 setNodes.insert(orderedNodes[i]);
2099 //MESSAGE("cellId = " << cellId);
2102 vtkIdType *nodes; // will refer to the point id's of the volume
2103 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
2106 //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};
2107 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,
2108 4, 0, 1, 5,16, 8,17,12, 5, 1, 2, 6,17, 9,18,13, 6, 2, 3, 7,18,10,19,14};
2109 for (int k = 0; k < 6; k++)
2112 for (int i = 0; i < 8; i++)
2113 tofind.insert(nodes[ids[8 * k + i]]);
2114 if (setNodes == tofind)
2116 for (int i = 0; i < 8; i++)
2117 orderedNodes[i] = nodes[ids[8 * k + i]];
2121 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
2122 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2123 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2126 void SMDS_DownQuadHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
2128 //ASSERT((cellId >=0)&& (cellId < _maxId));
2129 int *faces = &_cellIds[_nbDownCells * cellId];
2130 for (int i = 0; i < _nbDownCells; i++)
2134 faces[i] = lowCellId;
2137 if (faces[i] == lowCellId)
2143 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2144 * The ordering of the twenty points defining the quadratic hexahedron cell is point id's (0-7,8-19)
2145 * where point id's 0-7 are the eight corner vertices of the cube, followed by twelve mid-edge nodes (8-19).
2146 * Note that these mid-edge nodes lie on the edges defined by
2147 * (0,1), (1,2), (2,3), (3,0), (4,5), (5,6), (6,7), (7,4), (0,4), (1,5), (2,6), (3,7).
2148 * @see vtkQuadraticHexahedron.h in Filtering
2149 * @param cellId volumeId in vtkUnstructuredGrid
2150 * @param facesWithNodes vector of face descriptors to be filled
2152 void SMDS_DownQuadHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2154 // --- find point id's of the volume
2157 vtkIdType *nodes; // will refer to the point id's of the volume
2158 _grid->GetCellPoints(cellId, npts, nodes);
2160 // --- create all the ordered list of node id's for each face
2162 facesWithNodes.nbElems = 6;
2164 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2165 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2166 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2167 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2168 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
2169 facesWithNodes.elems[0].nodeIds[5] = nodes[9];
2170 facesWithNodes.elems[0].nodeIds[6] = nodes[10];
2171 facesWithNodes.elems[0].nodeIds[7] = nodes[11];
2172 facesWithNodes.elems[0].nbNodes = 8;
2173 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
2175 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2176 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2177 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2178 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2179 facesWithNodes.elems[1].nodeIds[4] = nodes[12];
2180 facesWithNodes.elems[1].nodeIds[5] = nodes[13];
2181 facesWithNodes.elems[1].nodeIds[6] = nodes[14];
2182 facesWithNodes.elems[1].nodeIds[7] = nodes[15];
2183 facesWithNodes.elems[1].nbNodes = 8;
2184 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
2186 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2187 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2188 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2189 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2190 facesWithNodes.elems[2].nodeIds[4] = nodes[8];
2191 facesWithNodes.elems[2].nodeIds[5] = nodes[17];
2192 facesWithNodes.elems[2].nodeIds[6] = nodes[12];
2193 facesWithNodes.elems[2].nodeIds[7] = nodes[16];
2194 facesWithNodes.elems[2].nbNodes = 8;
2195 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
2197 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2198 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2199 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2200 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2201 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
2202 facesWithNodes.elems[3].nodeIds[5] = nodes[18];
2203 facesWithNodes.elems[3].nodeIds[6] = nodes[13];
2204 facesWithNodes.elems[3].nodeIds[7] = nodes[17];
2205 facesWithNodes.elems[3].nbNodes = 8;
2206 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_QUAD;
2208 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2209 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2210 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2211 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2212 facesWithNodes.elems[4].nodeIds[4] = nodes[18];
2213 facesWithNodes.elems[4].nodeIds[5] = nodes[14];
2214 facesWithNodes.elems[4].nodeIds[6] = nodes[19];
2215 facesWithNodes.elems[4].nodeIds[7] = nodes[10];
2216 facesWithNodes.elems[4].nbNodes = 8;
2217 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_QUAD;
2219 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2220 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2221 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2222 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2223 facesWithNodes.elems[5].nodeIds[4] = nodes[19];
2224 facesWithNodes.elems[5].nodeIds[5] = nodes[15];
2225 facesWithNodes.elems[5].nodeIds[6] = nodes[16];
2226 facesWithNodes.elems[5].nodeIds[7] = nodes[11];
2227 facesWithNodes.elems[5].nbNodes = 8;
2228 facesWithNodes.elems[5].vtkType = VTK_QUADRATIC_QUAD;
2231 // ---------------------------------------------------------------------------