1 // Copyright (C) 2010-2013 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.
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>
35 // ---------------------------------------------------------------------------
37 vector<int> SMDS_Downward::_cellDimension;
39 /*! get the dimension of a cell (1,2,3 for 1D, 2D 3D) given the vtk cell type
41 * @param cellType vtk cell type @see vtkCellType.h
44 int SMDS_Downward::getCellDimension(unsigned char cellType)
46 if (_cellDimension.empty())
48 _cellDimension.resize(VTK_MAXTYPE + 1, 0);
49 _cellDimension[VTK_LINE] = 1;
50 _cellDimension[VTK_QUADRATIC_EDGE] = 1;
51 _cellDimension[VTK_TRIANGLE] = 2;
52 _cellDimension[VTK_QUADRATIC_TRIANGLE] = 2;
53 _cellDimension[VTK_QUAD] = 2;
54 _cellDimension[VTK_QUADRATIC_QUAD] = 2;
55 _cellDimension[VTK_BIQUADRATIC_QUAD] = 2;
56 _cellDimension[VTK_TETRA] = 3;
57 _cellDimension[VTK_QUADRATIC_TETRA] = 3;
58 _cellDimension[VTK_HEXAHEDRON] = 3;
59 _cellDimension[VTK_QUADRATIC_HEXAHEDRON] = 3;
60 _cellDimension[VTK_TRIQUADRATIC_HEXAHEDRON] = 3;
61 _cellDimension[VTK_WEDGE] = 3;
62 _cellDimension[VTK_QUADRATIC_WEDGE] = 3;
63 _cellDimension[VTK_PYRAMID] = 3;
64 _cellDimension[VTK_QUADRATIC_PYRAMID] = 3;
65 _cellDimension[VTK_HEXAGONAL_PRISM] = 3;
67 return _cellDimension[cellType];
70 // ---------------------------------------------------------------------------
72 /*! Generic constructor for all the downward connectivity structures (one per vtk cell type).
73 * The static structure for cell dimension is set only once.
74 * @param grid unstructured grid associated to the mesh.
75 * @param nbDownCells number of downward entities associated to this vtk type of cell.
78 SMDS_Downward::SMDS_Downward(SMDS_UnstructuredGrid *grid, int nbDownCells) :
79 _grid(grid), _nbDownCells(nbDownCells)
82 this->_cellIds.clear();
83 this->_cellTypes.clear();
84 if (_cellDimension.empty())
85 getCellDimension( VTK_LINE );
88 SMDS_Downward::~SMDS_Downward()
92 /*! Give or create an entry for downward connectivity structure relative to a cell.
93 * If the entry already exists, just return its id, otherwise, create it.
94 * The internal storage memory is allocated if needed.
95 * The SMDS_UnstructuredGrid::_cellIdToDownId vector is completed for vtkUnstructuredGrid cells.
96 * @param vtkId for a vtkUnstructuredGrid cell or -1 (default) for a created downward cell.
97 * @return the rank in downward[vtkType] structure.
99 int SMDS_Downward::addCell(int vtkId)
103 localId = _grid->CellIdToDownId(vtkId);
107 localId = this->_maxId;
109 this->allocate(_maxId);
112 this->_vtkCellIds[localId] = vtkId;
113 _grid->setCellIdToDownId(vtkId, localId);
115 this->initCell(localId);
119 /*! generic method do nothing. see derived methods
123 void SMDS_Downward::initCell(int cellId)
127 /*! Get the number of downward entities associated to a cell (always the same for a given vtk type of cell)
129 * @param cellId not used here.
132 int SMDS_Downward::getNumberOfDownCells(int cellId)
137 /*! get a pointer on the downward entities id's associated to a cell.
138 * @see SMDS_Downward::getNumberOfDownCells for the number of downward entities.
139 * @see SMDS_Downward::getDownTypes for the vtk cell types associated to the downward entities.
140 * @param cellId index of the cell in the downward structure relative to a given vtk cell type.
141 * @return table of downward entities id's.
143 const int* SMDS_Downward::getDownCells(int cellId)
145 //ASSERT((cellId >=0) && (cellId < _maxId));
146 return &_cellIds[_nbDownCells * cellId];
149 /*! get a list of vtk cell types associated to downward entities of a given cell, in the same order
150 * than the downward entities id's list (@see SMDS_Downward::getDownCells).
152 * @param cellId index of the cell in the downward structure relative to a vtk cell type.
153 * @return table of downward entities types.
155 const unsigned char* SMDS_Downward::getDownTypes(int cellId)
157 return &_cellTypes[0];
160 /*! add a downward entity of dimension n-1 (cell or node) to a given cell.
161 * Actual implementation is done in derived methods.
162 * @param cellId index of the parent cell (dimension n) in the downward structure relative to a vtk cell type.
163 * @param lowCellId index of the children cell to add (dimension n-1)
164 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
166 void SMDS_Downward::addDownCell(int cellId, int lowCellId, unsigned char aType)
168 ASSERT(0); // must be re-implemented in derived class
171 /*! add a downward entity of dimension n+1 to a given cell.
172 * Actual implementation is done in derived methods.
173 * @param cellId index of the children cell (dimension n) in the downward structure relative to a vtk cell type.
174 * @param upCellId index of the parent cell to add (dimension n+1)
175 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
177 void SMDS_Downward::addUpCell(int cellId, int upCellId, unsigned char aType)
179 ASSERT(0); // must be re-implemented in derived class
182 int SMDS_Downward::getNodeSet(int cellId, int* nodeSet)
187 // ---------------------------------------------------------------------------
189 SMDS_Down1D::SMDS_Down1D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
190 SMDS_Downward(grid, nbDownCells)
192 _upCellIdsVector.clear();
193 _upCellTypesVector.clear();
195 _upCellTypes.clear();
196 _upCellIndex.clear();
199 SMDS_Down1D::~SMDS_Down1D()
203 /*! clear vectors used to reference 2D cells containing the edge
207 void SMDS_Down1D::initCell(int cellId)
209 _upCellIdsVector[cellId].clear();
210 _upCellTypesVector[cellId].clear();
213 /*! Resize the downward connectivity storage vector if needed.
215 * @param nbElems total number of elements of the same type required
217 void SMDS_Down1D::allocate(int nbElems)
219 if (nbElems >= _vtkCellIds.size())
221 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
222 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
223 _upCellIdsVector.resize(nbElems + SMDS_Mesh::chunkSize);
224 _upCellTypesVector.resize(nbElems + SMDS_Mesh::chunkSize);
228 void SMDS_Down1D::compactStorage()
230 _cellIds.resize(_nbDownCells * _maxId);
231 _vtkCellIds.resize(_maxId);
234 for (int i = 0; i < _maxId; i++)
235 sizeUpCells += _upCellIdsVector[i].size();
236 _upCellIds.resize(sizeUpCells, -1);
237 _upCellTypes.resize(sizeUpCells);
238 _upCellIndex.resize(_maxId + 1, -1); // id and types of rank i correspond to [ _upCellIndex[i], _upCellIndex[i+1] [
240 for (int i = 0; i < _maxId; i++)
242 _upCellIndex[i] = current;
243 for (int j = 0; j < _upCellIdsVector[i].size(); j++)
245 _upCellIds[current] = _upCellIdsVector[i][j];
246 _upCellTypes[current] = _upCellTypesVector[i][j];
250 _upCellIndex[_maxId] = current;
252 _upCellIdsVector.clear();
253 _upCellTypesVector.clear();
256 void SMDS_Down1D::addUpCell(int cellId, int upCellId, unsigned char aType)
258 //ASSERT((cellId >=0) && (cellId < _maxId));
259 int nbFaces = _upCellIdsVector[cellId].size();
260 for (int i = 0; i < nbFaces; i++)
262 if ((_upCellIdsVector[cellId][i] == upCellId) && (_upCellTypesVector[cellId][i] == aType))
264 return; // already done
267 _upCellIdsVector[cellId].push_back(upCellId);
268 _upCellTypesVector[cellId].push_back(aType);
271 int SMDS_Down1D::getNumberOfUpCells(int cellId)
273 //ASSERT((cellId >=0) && (cellId < _maxId));
274 return _upCellIndex[cellId + 1] - _upCellIndex[cellId];
277 const int* SMDS_Down1D::getUpCells(int cellId)
279 //ASSERT((cellId >=0) && (cellId < _maxId));
280 return &_upCellIds[_upCellIndex[cellId]];
283 const unsigned char* SMDS_Down1D::getUpTypes(int cellId)
285 //ASSERT((cellId >=0) && (cellId < _maxId));
286 return &_upCellTypes[_upCellIndex[cellId]];
289 void SMDS_Down1D::getNodeIds(int cellId, std::set<int>& nodeSet)
291 for (int i = 0; i < _nbDownCells; i++)
292 nodeSet.insert(_cellIds[_nbDownCells * cellId + i]);
295 int SMDS_Down1D::getNodeSet(int cellId, int* nodeSet)
297 for (int i = 0; i < _nbDownCells; i++)
298 nodeSet[i] = _cellIds[_nbDownCells * cellId + i];
302 void SMDS_Down1D::setNodes(int cellId, int vtkId)
305 vtkIdType *pts; // will refer to the point id's of the face
306 _grid->GetCellPoints(vtkId, npts, pts);
307 // MESSAGE(vtkId << " " << npts << " " << _nbDownCells);
308 //ASSERT(npts == _nbDownCells);
309 for (int i = 0; i < npts; i++)
311 _cellIds[_nbDownCells * cellId + i] = pts[i];
315 void SMDS_Down1D::setNodes(int cellId, const int* nodeIds)
317 //ASSERT(nodeIds.size() == _nbDownCells);
318 for (int i = 0; i < _nbDownCells; i++)
320 _cellIds[_nbDownCells * cellId + i] = nodeIds[i];
324 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
325 * We keep in the list the cells that contains all the nodes, we keep only volumes and faces.
326 * @param cellId id of the edge in the downward structure
327 * @param vtkIds vector of vtk id's
328 * @return number of vtk cells (size of vector)
330 int SMDS_Down1D::computeVtkCells(int cellId, std::vector<int>& vtkIds)
334 // --- find all the cells the points belong to, and how many of the points belong to a given cell
336 int *pts = &_cellIds[_nbDownCells * cellId];
337 int ncells = this->computeVtkCells(pts, vtkIds);
341 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
343 * @param pts list of points id's defining an edge
344 * @param vtkIds vector of vtk id's
345 * @return number of vtk cells (size of vector)
347 int SMDS_Down1D::computeVtkCells(int *pts, std::vector<int>& vtkIds)
350 // --- find all the cells the points belong to, and how many of the points belong to a given cell
355 for (int i = 0; i < _nbDownCells; i++)
357 vtkIdType point = pts[i];
358 int numCells = _grid->GetLinks()->GetNcells(point);
359 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
360 for (int j = 0; j < numCells; j++)
362 int vtkCellId = cells[j];
364 for (int k = 0; k < cnt; k++)
366 if (cellIds[k] == vtkCellId)
375 cellIds[cnt] = vtkCellId;
377 // TODO ASSERT(cnt<1000);
383 // --- find the face and volume cells: they contains all the points and are of type volume or face
386 for (int i = 0; i < cnt; i++)
388 if (cellCnt[i] == _nbDownCells)
390 int vtkElemId = cellIds[i];
391 int vtkType = _grid->GetCellType(vtkElemId);
392 if (SMDS_Downward::getCellDimension(vtkType) > 1)
394 vtkIds.push_back(vtkElemId);
403 /*! Build the list of downward faces from a list of vtk cells.
405 * @param cellId id of the edge in the downward structure
406 * @param vtkIds vector of vtk id's
407 * @param downFaces vector of face id's in downward structures
408 * @param downTypes vector of face types
409 * @return number of downward faces
411 int SMDS_Down1D::computeFaces(int cellId, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
413 int *pts = &_cellIds[_nbDownCells * cellId];
414 int nbFaces = this->computeFaces(pts, vtkIds, nbcells, downFaces, downTypes);
418 /*! Build the list of downward faces from a list of vtk cells.
420 * @param pts list of points id's defining an edge
421 * @param vtkIds vector of vtk id's
422 * @param downFaces vector of face id's in downward structures
423 * @param downTypes vector of face types
424 * @return number of downward faces
426 int SMDS_Down1D::computeFaces(int* pts, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
429 for (int i = 0; i < nbcells; i++)
431 int vtkId = vtkIds[i];
432 int vtkType = _grid->GetCellType(vtkId);
433 if (SMDS_Downward::getCellDimension(vtkType) == 2)
435 int faceId = _grid->CellIdToDownId(vtkId);
436 downFaces[cnt] = faceId;
437 downTypes[cnt] = vtkType;
440 else if (SMDS_Downward::getCellDimension(vtkType) == 3)
442 int volId = _grid->CellIdToDownId(vtkId);
443 SMDS_Downward * downvol = _grid->getDownArray(vtkType);
444 //const int *downIds = downvol->getDownCells(volId);
445 const unsigned char* downTypesVol = downvol->getDownTypes(volId);
446 int nbFaces = downvol->getNumberOfDownCells(volId);
447 const int* faceIds = downvol->getDownCells(volId);
448 for (int n = 0; n < nbFaces; n++)
450 SMDS_Down2D *downFace = static_cast<SMDS_Down2D*> (_grid->getDownArray(downTypesVol[n]));
451 bool isInFace = downFace->isInFace(faceIds[n], pts, _nbDownCells);
454 bool alreadySet = false;
455 for (int k = 0; k < cnt; k++)
456 if (faceIds[n] == downFaces[k])
463 downFaces[cnt] = faceIds[n];
464 downTypes[cnt] = downTypesVol[n];
474 // ---------------------------------------------------------------------------
476 SMDS_Down2D::SMDS_Down2D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
477 SMDS_Downward(grid, nbDownCells)
480 _upCellTypes.clear();
485 SMDS_Down2D::~SMDS_Down2D()
489 int SMDS_Down2D::getNumberOfUpCells(int cellId)
492 if (_upCellIds[2 * cellId] >= 0)
494 if (_upCellIds[2 * cellId + 1] >= 0)
499 const int* SMDS_Down2D::getUpCells(int cellId)
501 //ASSERT((cellId >=0) && (cellId < _maxId));
502 return &_upCellIds[2 * cellId];
505 const unsigned char* SMDS_Down2D::getUpTypes(int cellId)
507 //ASSERT((cellId >=0) && (cellId < _maxId));
508 return &_upCellTypes[2 * cellId];
511 void SMDS_Down2D::getNodeIds(int cellId, std::set<int>& nodeSet)
513 for (int i = 0; i < _nbDownCells; i++)
515 int downCellId = _cellIds[_nbDownCells * cellId + i];
516 unsigned char cellType = _cellTypes[i];
517 this->_grid->getDownArray(cellType)->getNodeIds(downCellId, nodeSet);
521 /*! Find in vtkUnstructuredGrid the volumes containing a face already stored in vtkUnstructuredGrid.
522 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
523 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
524 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
525 * @param cellId the face cell id in vkUnstructuredGrid
526 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
527 * @return number of volumes (0, 1 or 2)
529 int SMDS_Down2D::computeVolumeIds(int cellId, int* ids)
531 // --- find point id's of the face
534 vtkIdType *pts; // will refer to the point id's of the face
535 _grid->GetCellPoints(cellId, npts, pts);
537 for (int i = 0; i < npts; i++)
538 nodes.push_back(pts[i]);
539 int nvol = this->computeVolumeIdsFromNodesFace(&nodes[0], npts, ids);
543 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
544 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
545 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
546 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
548 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
549 * @return number of volumes (0, 1 or 2)
551 int SMDS_Down2D::computeVolumeIds(ElemByNodesType& faceByNodes, int* ids)
553 int nvol = this->computeVolumeIdsFromNodesFace(&faceByNodes.nodeIds[0], faceByNodes.nbNodes, ids);
557 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
558 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
559 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
560 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
561 * @param pts array of vtk node id's
562 * @param npts number of nodes
564 * @return number of volumes (0, 1 or 2)
566 int SMDS_Down2D::computeVolumeIdsFromNodesFace(int* pts, int npts, int* ids)
569 // --- find all the cells the points belong to, and how many of the points belong to a given cell
574 for (int i = 0; i < npts; i++)
576 vtkIdType point = pts[i];
577 int numCells = _grid->GetLinks()->GetNcells(point);
578 //MESSAGE("cells pour " << i << " " << numCells);
579 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
580 for (int j = 0; j < numCells; j++)
582 int vtkCellId = cells[j];
584 for (int k = 0; k < cnt; k++)
586 if (cellIds[k] == vtkCellId)
595 cellIds[cnt] = vtkCellId;
597 // TODO ASSERT(cnt<1000);
603 // --- find the volume cells: they contains all the points and are of type volume
606 for (int i = 0; i < cnt; i++)
608 //MESSAGE("cell " << cellIds[i] << " points " << cellCnt[i]);
609 if (cellCnt[i] == npts)
611 int vtkElemId = cellIds[i];
612 int vtkType = _grid->GetCellType(vtkElemId);
613 if (SMDS_Downward::getCellDimension(vtkType) == 3)
615 ids[nvol] = vtkElemId; // store the volume id in given vector
626 void SMDS_Down2D::setTempNodes(int cellId, int vtkId)
629 vtkIdType *pts; // will refer to the point id's of the face
630 _grid->GetCellPoints(vtkId, npts, pts);
631 // MESSAGE(vtkId << " " << npts << " " << _nbNodes);
632 //ASSERT(npts == _nbNodes);
633 for (int i = 0; i < npts; i++)
635 _tempNodes[_nbNodes * cellId + i] = pts[i];
639 void SMDS_Down2D::setTempNodes(int cellId, ElemByNodesType& faceByNodes)
641 for (int i = 0; i < faceByNodes.nbNodes; i++)
642 _tempNodes[_nbNodes * cellId + i] = faceByNodes.nodeIds[i];
645 /*! Find if all the nodes belongs to the face.
647 * @param cellId the face cell Id
648 * @param nodeSet set of node id's to be found in the face list of nodes
651 bool SMDS_Down2D::isInFace(int cellId, int *pts, int npts)
654 int *nodes = &_tempNodes[_nbNodes * cellId];
655 for (int j = 0; j < npts; j++)
658 for (int i = 0; i < _nbNodes; i++)
660 if (nodes[i] == point)
667 return (nbFound == npts);
670 /*! Resize the downward connectivity storage vector if needed.
672 * @param nbElems total number of elements of the same type required
674 void SMDS_Down2D::allocate(int nbElems)
676 if (nbElems >= _vtkCellIds.size())
678 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
679 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
680 _upCellIds.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
681 _upCellTypes.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
682 _tempNodes.resize(_nbNodes * (nbElems + SMDS_Mesh::chunkSize), -1);
686 void SMDS_Down2D::compactStorage()
688 _cellIds.resize(_nbDownCells * _maxId);
689 _upCellIds.resize(2 * _maxId);
690 _upCellTypes.resize(2 * _maxId);
691 _vtkCellIds.resize(_maxId);
695 void SMDS_Down2D::addUpCell(int cellId, int upCellId, unsigned char aType)
697 //ASSERT((cellId >=0)&& (cellId < _maxId));
698 int *vols = &_upCellIds[2 * cellId];
699 unsigned char *types = &_upCellTypes[2 * cellId];
700 for (int i = 0; i < 2; i++)
704 vols[i] = upCellId; // use non affected volume
708 if ((vols[i] == upCellId) && (types[i] == aType)) // already done
714 int SMDS_Down2D::getNodeSet(int cellId, int* nodeSet)
716 for (int i = 0; i < _nbNodes; i++)
717 nodeSet[i] = _tempNodes[_nbNodes * cellId + i];
721 int SMDS_Down2D::FindEdgeByNodes(int cellId, ElemByNodesType& edgeByNodes)
723 int *edges = &_cellIds[_nbDownCells * cellId];
724 for (int i = 0; i < _nbDownCells; i++)
726 if ((edges[i] >= 0) && (edgeByNodes.vtkType == _cellTypes[i]))
729 int npts = this->_grid->getDownArray(edgeByNodes.vtkType)->getNodeSet(edges[i], nodeSet);
731 for (int j = 0; j < npts; j++)
733 int point = edgeByNodes.nodeIds[j];
735 for (int k = 0; k < npts; k++)
737 if (nodeSet[k] == point)
753 // ---------------------------------------------------------------------------
755 SMDS_Down3D::SMDS_Down3D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
756 SMDS_Downward(grid, nbDownCells)
760 SMDS_Down3D::~SMDS_Down3D()
764 void SMDS_Down3D::allocate(int nbElems)
766 if (nbElems >= _vtkCellIds.size())
768 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
769 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
773 void SMDS_Down3D::compactStorage()
775 // nothing to do, size was known before
778 int SMDS_Down3D::getNumberOfUpCells(int cellId)
783 const int* SMDS_Down3D::getUpCells(int cellId)
788 const unsigned char* SMDS_Down3D::getUpTypes(int cellId)
793 void SMDS_Down3D::getNodeIds(int cellId, std::set<int>& nodeSet)
795 int vtkId = this->_vtkCellIds[cellId];
797 vtkIdType *nodes; // will refer to the point id's of the volume
798 _grid->GetCellPoints(vtkId, npts, nodes);
799 for (int i = 0; i < npts; i++)
800 nodeSet.insert(nodes[i]);
803 int SMDS_Down3D::FindFaceByNodes(int cellId, ElemByNodesType& faceByNodes)
805 int *faces = &_cellIds[_nbDownCells * cellId];
808 for (int i = 0; i < _nbDownCells; i++)
810 if ((faces[i] >= 0) && (faceByNodes.vtkType == _cellTypes[i]))
813 npoints = faceByNodes.nbNodes;
816 int npts = this->_grid->getDownArray(faceByNodes.vtkType)->getNodeSet(faces[i], nodeSet);
818 continue; // skip this face
820 for (int j = 0; j < npts; j++)
822 int point = faceByNodes.nodeIds[j];
824 for (int k = 0; k < npts; k++)
826 if (nodeSet[k] == point)
829 break; // point j is in the 2 faces, skip remaining k values
833 break; // point j is not in the 2 faces, skip the remaining tests
842 // ---------------------------------------------------------------------------
844 SMDS_DownEdge::SMDS_DownEdge(SMDS_UnstructuredGrid *grid) :
847 _cellTypes.push_back(VTK_VERTEX);
848 _cellTypes.push_back(VTK_VERTEX);
851 SMDS_DownEdge::~SMDS_DownEdge()
855 // ---------------------------------------------------------------------------
857 SMDS_DownQuadEdge::SMDS_DownQuadEdge(SMDS_UnstructuredGrid *grid) :
860 _cellTypes.push_back(VTK_VERTEX);
861 _cellTypes.push_back(VTK_VERTEX);
862 _cellTypes.push_back(VTK_VERTEX);
865 SMDS_DownQuadEdge::~SMDS_DownQuadEdge()
869 // ---------------------------------------------------------------------------
871 SMDS_DownTriangle::SMDS_DownTriangle(SMDS_UnstructuredGrid *grid) :
874 _cellTypes.push_back(VTK_LINE);
875 _cellTypes.push_back(VTK_LINE);
876 _cellTypes.push_back(VTK_LINE);
880 SMDS_DownTriangle::~SMDS_DownTriangle()
884 void SMDS_DownTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
886 int *nodes = &_tempNodes[_nbNodes * cellId];
887 edgesWithNodes.nbElems = 3;
889 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
890 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
891 edgesWithNodes.elems[0].nbNodes = 2;
892 edgesWithNodes.elems[0].vtkType = VTK_LINE;
894 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
895 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
896 edgesWithNodes.elems[1].nbNodes = 2;
897 edgesWithNodes.elems[1].vtkType = VTK_LINE;
899 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
900 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
901 edgesWithNodes.elems[2].nbNodes = 2;
902 edgesWithNodes.elems[2].vtkType = VTK_LINE;
905 void SMDS_DownTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
907 //ASSERT((cellId >=0)&& (cellId < _maxId));
908 //ASSERT(aType == VTK_LINE);
909 int *faces = &_cellIds[_nbDownCells * cellId];
910 for (int i = 0; i < _nbDownCells; i++)
914 faces[i] = lowCellId;
917 if (faces[i] == lowCellId)
923 // ---------------------------------------------------------------------------
925 SMDS_DownQuadTriangle::SMDS_DownQuadTriangle(SMDS_UnstructuredGrid *grid) :
928 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
929 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
930 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
934 SMDS_DownQuadTriangle::~SMDS_DownQuadTriangle()
938 void SMDS_DownQuadTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
940 int *nodes = &_tempNodes[_nbNodes * cellId];
941 edgesWithNodes.nbElems = 3;
943 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
944 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
945 edgesWithNodes.elems[0].nodeIds[2] = nodes[3];
946 edgesWithNodes.elems[0].nbNodes = 3;
947 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
949 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
950 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
951 edgesWithNodes.elems[1].nodeIds[2] = nodes[4];
952 edgesWithNodes.elems[1].nbNodes = 3;
953 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
955 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
956 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
957 edgesWithNodes.elems[2].nodeIds[2] = nodes[5];
958 edgesWithNodes.elems[2].nbNodes = 3;
959 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
962 void SMDS_DownQuadTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
964 //ASSERT((cellId >=0)&& (cellId < _maxId));
965 //ASSERT(aType == VTK_QUADRATIC_EDGE);
966 int *edges = &_cellIds[_nbDownCells * cellId];
967 for (int i = 0; i < _nbDownCells; i++)
971 edges[i] = lowCellId;
974 if (edges[i] == lowCellId)
980 // ---------------------------------------------------------------------------
982 SMDS_DownQuadrangle::SMDS_DownQuadrangle(SMDS_UnstructuredGrid *grid) :
985 _cellTypes.push_back(VTK_LINE);
986 _cellTypes.push_back(VTK_LINE);
987 _cellTypes.push_back(VTK_LINE);
988 _cellTypes.push_back(VTK_LINE);
992 SMDS_DownQuadrangle::~SMDS_DownQuadrangle()
996 void SMDS_DownQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
998 int *nodes = &_tempNodes[_nbNodes * cellId];
999 edgesWithNodes.nbElems = 4;
1001 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1002 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1003 edgesWithNodes.elems[0].nbNodes = 2;
1004 edgesWithNodes.elems[0].vtkType = VTK_LINE;
1006 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1007 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1008 edgesWithNodes.elems[1].nbNodes = 2;
1009 edgesWithNodes.elems[1].vtkType = VTK_LINE;
1011 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1012 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1013 edgesWithNodes.elems[2].nbNodes = 2;
1014 edgesWithNodes.elems[2].vtkType = VTK_LINE;
1016 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1017 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1018 edgesWithNodes.elems[3].nbNodes = 2;
1019 edgesWithNodes.elems[3].vtkType = VTK_LINE;
1022 void SMDS_DownQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1024 //ASSERT((cellId >=0)&& (cellId < _maxId));
1025 //ASSERT(aType == VTK_LINE);
1026 int *faces = &_cellIds[_nbDownCells * cellId];
1027 for (int i = 0; i < _nbDownCells; i++)
1031 faces[i] = lowCellId;
1034 if (faces[i] == lowCellId)
1040 // ---------------------------------------------------------------------------
1042 SMDS_DownQuadQuadrangle::SMDS_DownQuadQuadrangle(SMDS_UnstructuredGrid *grid) :
1043 SMDS_Down2D(grid, 4)
1045 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1046 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1047 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1048 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1052 SMDS_DownQuadQuadrangle::~SMDS_DownQuadQuadrangle()
1056 void SMDS_DownQuadQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
1058 int *nodes = &_tempNodes[_nbNodes * cellId];
1059 edgesWithNodes.nbElems = 4;
1061 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1062 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1063 edgesWithNodes.elems[0].nodeIds[2] = nodes[4];
1064 edgesWithNodes.elems[0].nbNodes = 3;
1065 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
1067 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1068 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1069 edgesWithNodes.elems[1].nodeIds[2] = nodes[5];
1070 edgesWithNodes.elems[1].nbNodes = 3;
1071 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
1073 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1074 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1075 edgesWithNodes.elems[2].nodeIds[2] = nodes[6];
1076 edgesWithNodes.elems[2].nbNodes = 3;
1077 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
1079 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1080 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1081 edgesWithNodes.elems[3].nodeIds[2] = nodes[7];
1082 edgesWithNodes.elems[3].nbNodes = 3;
1083 edgesWithNodes.elems[3].vtkType = VTK_QUADRATIC_EDGE;
1086 void SMDS_DownQuadQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1088 //ASSERT((cellId >=0)&& (cellId < _maxId));
1089 //ASSERT(aType == VTK_QUADRATIC_EDGE);
1090 int *faces = &_cellIds[_nbDownCells * cellId];
1091 for (int i = 0; i < _nbDownCells; i++)
1095 faces[i] = lowCellId;
1098 if (faces[i] == lowCellId)
1104 // ---------------------------------------------------------------------------
1106 SMDS_DownTetra::SMDS_DownTetra(SMDS_UnstructuredGrid *grid) :
1107 SMDS_Down3D(grid, 4)
1109 _cellTypes.push_back(VTK_TRIANGLE);
1110 _cellTypes.push_back(VTK_TRIANGLE);
1111 _cellTypes.push_back(VTK_TRIANGLE);
1112 _cellTypes.push_back(VTK_TRIANGLE);
1115 SMDS_DownTetra::~SMDS_DownTetra()
1119 void SMDS_DownTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1123 for (int i = 0; i < orderedNodes.size(); i++)
1124 setNodes.insert(orderedNodes[i]);
1125 //MESSAGE("cellId = " << cellId);
1128 vtkIdType *nodes; // will refer to the point id's of the volume
1129 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1132 int ids[12] = { 0, 1, 2, 0, 3, 1, 2, 3, 0, 1, 3, 2 };
1133 //int ids[12] = { 2, 1, 0, 1, 3, 0, 0, 3, 2, 2, 3, 1 };
1134 for (int k = 0; k < 4; k++)
1137 for (int i = 0; i < 3; i++)
1138 tofind.insert(nodes[ids[3 * k + i]]);
1139 if (setNodes == tofind)
1141 for (int i = 0; i < 3; i++)
1142 orderedNodes[i] = nodes[ids[3 * k + i]];
1146 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1147 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1148 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1151 void SMDS_DownTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1153 //ASSERT((cellId >=0)&& (cellId < _maxId));
1154 //ASSERT(aType == VTK_TRIANGLE);
1155 int *faces = &_cellIds[_nbDownCells * cellId];
1156 for (int i = 0; i < _nbDownCells; i++)
1160 faces[i] = lowCellId;
1163 if (faces[i] == lowCellId)
1169 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1170 * The linear tetrahedron is defined by four points.
1171 * @see vtkTetra.h in Filtering.
1172 * @param cellId volumeId in vtkUnstructuredGrid
1173 * @param facesWithNodes vector of face descriptors to be filled
1175 void SMDS_DownTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1177 // --- find point id's of the volume
1180 vtkIdType *nodes; // will refer to the point id's of the volume
1181 _grid->GetCellPoints(cellId, npts, nodes);
1183 // --- create all the ordered list of node id's for each face
1185 facesWithNodes.nbElems = 4;
1187 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1188 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1189 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1190 facesWithNodes.elems[0].nbNodes = 3;
1191 facesWithNodes.elems[0].vtkType = VTK_TRIANGLE;
1193 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1194 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1195 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1196 facesWithNodes.elems[1].nbNodes = 3;
1197 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1199 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1200 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1201 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1202 facesWithNodes.elems[2].nbNodes = 3;
1203 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1205 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1206 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1207 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1208 facesWithNodes.elems[3].nbNodes = 3;
1209 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1212 // ---------------------------------------------------------------------------
1214 SMDS_DownQuadTetra::SMDS_DownQuadTetra(SMDS_UnstructuredGrid *grid) :
1215 SMDS_Down3D(grid, 4)
1217 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1218 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1219 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1220 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1223 SMDS_DownQuadTetra::~SMDS_DownQuadTetra()
1227 void SMDS_DownQuadTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1231 for (int i = 0; i < orderedNodes.size(); i++)
1232 setNodes.insert(orderedNodes[i]);
1233 //MESSAGE("cellId = " << cellId);
1236 vtkIdType *nodes; // will refer to the point id's of the volume
1237 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1240 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 };
1241 //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 };
1242 for (int k = 0; k < 4; k++)
1245 for (int i = 0; i < 6; i++)
1246 tofind.insert(nodes[ids[6 * k + i]]);
1247 if (setNodes == tofind)
1249 for (int i = 0; i < 6; i++)
1250 orderedNodes[i] = nodes[ids[6 * k + i]];
1254 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1255 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1256 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1259 void SMDS_DownQuadTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1261 //ASSERT((cellId >=0)&& (cellId < _maxId));
1262 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1263 int *faces = &_cellIds[_nbDownCells * cellId];
1264 for (int i = 0; i < _nbDownCells; i++)
1268 faces[i] = lowCellId;
1271 if (faces[i] == lowCellId)
1277 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1278 * The ordering of the ten points defining the quadratic tetrahedron cell is point id's (0-3,4-9)
1279 * where id's 0-3 are the four tetrahedron vertices;
1280 * 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).
1281 * @see vtkQuadraticTetra.h in Filtering.
1282 * @param cellId volumeId in vtkUnstructuredGrid
1283 * @param facesWithNodes vector of face descriptors to be filled
1285 void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1287 // --- find point id's of the volume
1290 vtkIdType *nodes; // will refer to the point id's of the volume
1291 _grid->GetCellPoints(cellId, npts, nodes);
1293 // --- create all the ordered list of node id's for each face
1295 facesWithNodes.nbElems = 4;
1297 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1298 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1299 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1300 facesWithNodes.elems[0].nodeIds[3] = nodes[4];
1301 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1302 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1303 facesWithNodes.elems[0].nbNodes = 6;
1304 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_TRIANGLE;
1306 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1307 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1308 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1309 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1310 facesWithNodes.elems[1].nodeIds[4] = nodes[8];
1311 facesWithNodes.elems[1].nodeIds[5] = nodes[7];
1312 facesWithNodes.elems[1].nbNodes = 6;
1313 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1315 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1316 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1317 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1318 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1319 facesWithNodes.elems[2].nodeIds[4] = nodes[9];
1320 facesWithNodes.elems[2].nodeIds[5] = nodes[7];
1321 facesWithNodes.elems[2].nbNodes = 6;
1322 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1324 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1325 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1326 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1327 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
1328 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
1329 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1330 facesWithNodes.elems[3].nbNodes = 6;
1331 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1334 // ---------------------------------------------------------------------------
1336 SMDS_DownPyramid::SMDS_DownPyramid(SMDS_UnstructuredGrid *grid) :
1337 SMDS_Down3D(grid, 5)
1339 _cellTypes.push_back(VTK_QUAD);
1340 _cellTypes.push_back(VTK_TRIANGLE);
1341 _cellTypes.push_back(VTK_TRIANGLE);
1342 _cellTypes.push_back(VTK_TRIANGLE);
1343 _cellTypes.push_back(VTK_TRIANGLE);
1346 SMDS_DownPyramid::~SMDS_DownPyramid()
1350 void SMDS_DownPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1354 for (int i = 0; i < orderedNodes.size(); i++)
1355 setNodes.insert(orderedNodes[i]);
1356 //MESSAGE("cellId = " << cellId);
1359 vtkIdType *nodes; // will refer to the point id's of the volume
1360 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1363 int ids[16] = { 0, 1, 2, 3, 0, 3, 4, 3, 2, 4, 2, 1, 4, 1, 0, 4 };
1365 // Quadrangular face
1367 for (int i = 0; i < 4; i++)
1368 tofind.insert(nodes[ids[i]]);
1369 if (setNodes == tofind)
1371 for (int i = 0; i < 4; i++)
1372 orderedNodes[i] = nodes[ids[i]];
1376 for (int k = 0; k < 4; k++)
1379 for (int i = 0; i < 3; i++)
1380 tofind.insert(nodes[ids[4 + 3 * k + i]]);
1381 if (setNodes == tofind)
1383 for (int i = 0; i < 3; i++)
1384 orderedNodes[i] = nodes[ids[4 + 3 * k + i]];
1388 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1389 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1390 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1393 void SMDS_DownPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1395 //ASSERT((cellId >=0) && (cellId < _maxId));
1396 int *faces = &_cellIds[_nbDownCells * cellId];
1397 if (aType == VTK_QUAD)
1401 faces[0] = lowCellId;
1404 if (faces[0] == lowCellId)
1409 //ASSERT(aType == VTK_TRIANGLE);
1410 for (int i = 1; i < _nbDownCells; i++)
1414 faces[i] = lowCellId;
1417 if (faces[i] == lowCellId)
1424 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1425 * The pyramid is defined by the five points (0-4) where (0,1,2,3) is the base of the pyramid which,
1426 * using the right hand rule, forms a quadrilateral whose normal points in the direction of the
1427 * pyramid apex at vertex #4.
1428 * @see vtkPyramid.h in Filtering.
1429 * @param cellId volumeId in vtkUnstructuredGrid
1430 * @param facesWithNodes vector of face descriptors to be filled
1432 void SMDS_DownPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1434 // --- find point id's of the volume
1437 vtkIdType *nodes; // will refer to the point id's of the volume
1438 _grid->GetCellPoints(cellId, npts, nodes);
1440 // --- create all the ordered list of node id's for each face
1442 facesWithNodes.nbElems = 5;
1444 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1445 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1446 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1447 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1448 facesWithNodes.elems[0].nbNodes = 4;
1449 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1451 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1452 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1453 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1454 facesWithNodes.elems[1].nbNodes = 3;
1455 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1457 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1458 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1459 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1460 facesWithNodes.elems[2].nbNodes = 3;
1461 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1463 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1464 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1465 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1466 facesWithNodes.elems[3].nbNodes = 3;
1467 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1469 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1470 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1471 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1472 facesWithNodes.elems[4].nbNodes = 3;
1473 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1476 // ---------------------------------------------------------------------------
1478 SMDS_DownQuadPyramid::SMDS_DownQuadPyramid(SMDS_UnstructuredGrid *grid) :
1479 SMDS_Down3D(grid, 5)
1481 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1482 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1483 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1484 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1485 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1488 SMDS_DownQuadPyramid::~SMDS_DownQuadPyramid()
1492 void SMDS_DownQuadPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1494 // MESSAGE("SMDS_DownQuadPyramid::getOrderedNodesOfFace cellId = " << cellId);
1497 for (int i = 0; i < orderedNodes.size(); i++)
1498 setNodes.insert(orderedNodes[i]);
1499 //MESSAGE("cellId = " << cellId);
1502 vtkIdType *nodes; // will refer to the point id's of the volume
1503 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1506 int ids[32] = { 0, 1, 2, 3, 5, 6, 7, 8,
1507 0, 3, 4, 8, 12, 9, 3, 2, 4, 7 , 11, 12, 2, 1, 4, 6, 10, 11, 1, 0, 4, 5, 9, 10 };
1509 // Quadrangular face
1511 for (int i = 0; i < 8; i++)
1512 tofind.insert(nodes[ids[i]]);
1513 if (setNodes == tofind)
1515 for (int i = 0; i < 8; i++)
1516 orderedNodes[i] = nodes[ids[i]];
1520 for (int k = 0; k < 4; k++)
1523 for (int i = 0; i < 6; i++)
1524 tofind.insert(nodes[ids[8 + 6 * k + i]]);
1525 if (setNodes == tofind)
1527 for (int i = 0; i < 6; i++)
1528 orderedNodes[i] = nodes[ids[8 + 6 * k + i]];
1532 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1533 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
1534 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1537 void SMDS_DownQuadPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1539 //ASSERT((cellId >=0) && (cellId < _maxId));
1540 int *faces = &_cellIds[_nbDownCells * cellId];
1541 if (aType == VTK_QUADRATIC_QUAD)
1545 faces[0] = lowCellId;
1548 if (faces[0] == lowCellId)
1553 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1554 for (int i = 1; i < _nbDownCells; i++)
1558 faces[i] = lowCellId;
1561 if (faces[i] == lowCellId)
1568 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1569 * The ordering of the thirteen points defining the quadratic pyramid cell is point id's (0-4,5-12)
1570 * where point id's 0-4 are the five corner vertices of the pyramid; followed
1571 * by eight mid-edge nodes (5-12). Note that these mid-edge nodes lie on the edges defined by
1572 * 5(0,1), 6(1,2), 7(2,3), 8(3,0), 9(0,4), 10(1,4), 11(2,4), 12(3,4).
1573 * @see vtkQuadraticPyramid.h in Filtering.
1574 * @param cellId volumeId in vtkUnstructuredGrid
1575 * @param facesWithNodes vector of face descriptors to be filled
1577 void SMDS_DownQuadPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1579 // --- find point id's of the volume
1582 vtkIdType *nodes; // will refer to the point id's of the volume
1583 _grid->GetCellPoints(cellId, npts, nodes);
1585 // --- create all the ordered list of node id's for each face
1587 facesWithNodes.nbElems = 5;
1589 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1590 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1591 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1592 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1593 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1594 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1595 facesWithNodes.elems[0].nodeIds[6] = nodes[7];
1596 facesWithNodes.elems[0].nodeIds[7] = nodes[8];
1597 facesWithNodes.elems[0].nbNodes = 8;
1598 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1600 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1601 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1602 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1603 facesWithNodes.elems[1].nodeIds[3] = nodes[5];
1604 facesWithNodes.elems[1].nodeIds[4] = nodes[10];
1605 facesWithNodes.elems[1].nodeIds[5] = nodes[9];
1606 facesWithNodes.elems[1].nbNodes = 6;
1607 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1609 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1610 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1611 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1612 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1613 facesWithNodes.elems[2].nodeIds[4] = nodes[11];
1614 facesWithNodes.elems[2].nodeIds[5] = nodes[10];
1615 facesWithNodes.elems[2].nbNodes = 6;
1616 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1618 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1619 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1620 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1621 facesWithNodes.elems[3].nodeIds[3] = nodes[7];
1622 facesWithNodes.elems[3].nodeIds[4] = nodes[12];
1623 facesWithNodes.elems[3].nodeIds[5] = nodes[11];
1624 facesWithNodes.elems[3].nbNodes = 6;
1625 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1627 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1628 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1629 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1630 facesWithNodes.elems[4].nodeIds[3] = nodes[8];
1631 facesWithNodes.elems[4].nodeIds[4] = nodes[9];
1632 facesWithNodes.elems[4].nodeIds[5] = nodes[12];
1633 facesWithNodes.elems[4].nbNodes = 6;
1634 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1637 // ---------------------------------------------------------------------------
1639 SMDS_DownPenta::SMDS_DownPenta(SMDS_UnstructuredGrid *grid) :
1640 SMDS_Down3D(grid, 5)
1642 _cellTypes.push_back(VTK_QUAD);
1643 _cellTypes.push_back(VTK_QUAD);
1644 _cellTypes.push_back(VTK_QUAD);
1645 _cellTypes.push_back(VTK_TRIANGLE);
1646 _cellTypes.push_back(VTK_TRIANGLE);
1649 SMDS_DownPenta::~SMDS_DownPenta()
1653 void SMDS_DownPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1657 for (int i = 0; i < orderedNodes.size(); i++)
1658 setNodes.insert(orderedNodes[i]);
1659 //MESSAGE("cellId = " << cellId);
1662 vtkIdType *nodes; // will refer to the point id's of the volume
1663 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1666 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1667 int ids[18] = { 0, 1, 2, 3, 5, 4, 0, 3, 4, 1, 1, 4, 5, 2, 2, 5, 3, 0 };
1670 for (int k = 0; k < 2; k++)
1673 for (int i = 0; i < 3; i++)
1674 tofind.insert(nodes[ids[3 * k + i]]);
1675 if (setNodes == tofind)
1677 for (int i = 0; i < 3; i++)
1678 orderedNodes[i] = nodes[ids[3 * k + i]];
1682 // Quadrangular faces
1683 for (int k = 0; k < 3; k++)
1686 for (int i = 0; i < 4; i++)
1687 tofind.insert(nodes[ids[6 + 4 * k + i]]);
1688 if (setNodes == tofind)
1690 for (int i = 0; i < 4; i++)
1691 orderedNodes[i] = nodes[ids[6 + 4 * k + i]];
1695 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1696 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1697 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1700 void SMDS_DownPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1702 //ASSERT((cellId >=0) && (cellId < _maxId));
1703 int *faces = &_cellIds[_nbDownCells * cellId];
1704 if (aType == VTK_QUAD)
1705 for (int i = 0; i < 3; i++)
1709 faces[i] = lowCellId;
1712 if (faces[i] == lowCellId)
1717 //ASSERT(aType == VTK_TRIANGLE);
1718 for (int i = 3; i < _nbDownCells; i++)
1722 faces[i] = lowCellId;
1725 if (faces[i] == lowCellId)
1732 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's.
1733 * A wedge or pentahedron consists of two triangular and three quadrilateral faces
1734 * and is defined by the six points (0-5) where (0,1,2) is the base of the wedge which,
1735 * using the right hand rule, forms a triangle whose normal points outward
1736 * (away from the triangular face (3,4,5)).
1737 * @see vtkWedge.h in Filtering
1738 * @param cellId volumeId in vtkUnstructuredGrid
1739 * @param facesWithNodes vector of face descriptors to be filled
1741 void SMDS_DownPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1743 // --- find point id's of the volume
1746 vtkIdType *nodes; // will refer to the point id's of the volume
1747 _grid->GetCellPoints(cellId, npts, nodes);
1749 // --- create all the ordered list of node id's for each face
1751 facesWithNodes.nbElems = 5;
1753 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1754 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1755 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1756 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1757 facesWithNodes.elems[0].nbNodes = 4;
1758 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1760 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1761 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1762 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1763 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1764 facesWithNodes.elems[1].nbNodes = 4;
1765 facesWithNodes.elems[1].vtkType = VTK_QUAD;
1767 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1768 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1769 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1770 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1771 facesWithNodes.elems[2].nbNodes = 4;
1772 facesWithNodes.elems[2].vtkType = VTK_QUAD;
1774 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1775 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1776 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1777 facesWithNodes.elems[3].nbNodes = 3;
1778 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1780 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1781 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1782 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1783 facesWithNodes.elems[4].nbNodes = 3;
1784 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1787 // ---------------------------------------------------------------------------
1789 SMDS_DownQuadPenta::SMDS_DownQuadPenta(SMDS_UnstructuredGrid *grid) :
1790 SMDS_Down3D(grid, 5)
1792 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1793 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1794 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1795 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1796 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1799 SMDS_DownQuadPenta::~SMDS_DownQuadPenta()
1803 void SMDS_DownQuadPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1807 for (int i = 0; i < orderedNodes.size(); i++)
1808 setNodes.insert(orderedNodes[i]);
1809 //MESSAGE("cellId = " << cellId);
1812 vtkIdType *nodes; // will refer to the point id's of the volume
1813 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1816 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1817 int ids[36] = { 0, 1, 2, 6, 7, 8, 3, 5, 4, 11, 10, 9,
1818 0, 3, 4, 1, 12, 9, 13, 6, 1, 4, 5, 2, 13, 10, 14, 7, 2, 5, 3, 0, 14, 11, 12, 8 };
1821 for (int k = 0; k < 2; k++)
1824 for (int i = 0; i < 6; i++)
1825 tofind.insert(nodes[ids[6 * k + i]]);
1826 if (setNodes == tofind)
1828 for (int i = 0; i < 6; i++)
1829 orderedNodes[i] = nodes[ids[6 * k + i]];
1833 // Quadrangular faces
1834 for (int k = 0; k < 3; k++)
1837 for (int i = 0; i < 8; i++)
1838 tofind.insert(nodes[ids[12 + 8 * k + i]]);
1839 if (setNodes == tofind)
1841 for (int i = 0; i < 8; i++)
1842 orderedNodes[i] = nodes[ids[12 + 8 * k + i]];
1846 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1847 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1848 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1851 void SMDS_DownQuadPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1853 //ASSERT((cellId >=0) && (cellId < _maxId));
1854 int *faces = &_cellIds[_nbDownCells * cellId];
1855 if (aType == VTK_QUADRATIC_QUAD)
1856 for (int i = 0; i < 3; i++)
1860 faces[i] = lowCellId;
1863 if (faces[i] == lowCellId)
1868 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1869 for (int i = 3; i < _nbDownCells; i++)
1873 faces[i] = lowCellId;
1876 if (faces[i] == lowCellId)
1883 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1884 * The quadratic wedge (or pentahedron) is defined by fifteen points.
1885 * The ordering of the fifteen points defining the cell is point id's (0-5,6-14)
1886 * where point id's 0-5 are the six corner vertices of the wedge, followed by
1887 * nine mid-edge nodes (6-14). Note that these mid-edge nodes lie on the edges defined by
1888 * (0,1), (1,2), (2,0), (3,4), (4,5), (5,3), (0,3), (1,4), (2,5).
1889 * @see vtkQuadraticWedge.h in Filtering
1890 * @param cellId volumeId in vtkUnstructuredGrid
1891 * @param facesWithNodes vector of face descriptors to be filled
1893 void SMDS_DownQuadPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1895 // --- find point id's of the volume
1898 vtkIdType *nodes; // will refer to the point id's of the volume
1899 _grid->GetCellPoints(cellId, npts, nodes);
1901 // --- create all the ordered list of node id's for each face
1903 facesWithNodes.nbElems = 5;
1905 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1906 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1907 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1908 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1909 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
1910 facesWithNodes.elems[0].nodeIds[5] = nodes[14];
1911 facesWithNodes.elems[0].nodeIds[6] = nodes[11];
1912 facesWithNodes.elems[0].nodeIds[7] = nodes[12];
1913 facesWithNodes.elems[0].nbNodes = 8;
1914 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1916 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1917 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1918 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1919 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1920 facesWithNodes.elems[1].nodeIds[4] = nodes[7];
1921 facesWithNodes.elems[1].nodeIds[5] = nodes[14];
1922 facesWithNodes.elems[1].nodeIds[6] = nodes[10];
1923 facesWithNodes.elems[1].nodeIds[7] = nodes[13];
1924 facesWithNodes.elems[1].nbNodes = 8;
1925 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
1927 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1928 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1929 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1930 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1931 facesWithNodes.elems[2].nodeIds[4] = nodes[6];
1932 facesWithNodes.elems[2].nodeIds[5] = nodes[13];
1933 facesWithNodes.elems[2].nodeIds[6] = nodes[9];
1934 facesWithNodes.elems[2].nodeIds[7] = nodes[12];
1935 facesWithNodes.elems[2].nbNodes = 8;
1936 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
1938 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1939 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1940 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1941 facesWithNodes.elems[3].nodeIds[3] = nodes[6];
1942 facesWithNodes.elems[3].nodeIds[4] = nodes[7];
1943 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1944 facesWithNodes.elems[3].nbNodes = 6;
1945 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1947 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1948 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1949 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1950 facesWithNodes.elems[4].nodeIds[3] = nodes[9];
1951 facesWithNodes.elems[4].nodeIds[4] = nodes[10];
1952 facesWithNodes.elems[4].nodeIds[5] = nodes[11];
1953 facesWithNodes.elems[4].nbNodes = 6;
1954 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1957 // ---------------------------------------------------------------------------
1959 SMDS_DownHexa::SMDS_DownHexa(SMDS_UnstructuredGrid *grid) :
1960 SMDS_Down3D(grid, 6)
1962 _cellTypes.push_back(VTK_QUAD);
1963 _cellTypes.push_back(VTK_QUAD);
1964 _cellTypes.push_back(VTK_QUAD);
1965 _cellTypes.push_back(VTK_QUAD);
1966 _cellTypes.push_back(VTK_QUAD);
1967 _cellTypes.push_back(VTK_QUAD);
1970 SMDS_DownHexa::~SMDS_DownHexa()
1974 void SMDS_DownHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1978 for (int i = 0; i < orderedNodes.size(); i++)
1979 setNodes.insert(orderedNodes[i]);
1980 //MESSAGE("cellId = " << cellId);
1983 vtkIdType *nodes; // will refer to the point id's of the volume
1984 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1987 //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};
1988 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};
1989 for (int k = 0; k < 6; k++) // loop on the 6 faces
1992 for (int i = 0; i < 4; i++)
1993 tofind.insert(nodes[ids[4 * k + i]]); // node ids of the face i
1994 if (setNodes == tofind)
1996 for (int i = 0; i < 4; i++)
1997 orderedNodes[i] = nodes[ids[4 * k + i]];
2001 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
2002 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2003 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2004 MESSAGE(nodes[4] << " " << nodes[5] << " " << nodes[6] << " " << nodes[7]);
2007 void SMDS_DownHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
2009 //ASSERT((cellId >=0)&& (cellId < _maxId));
2010 int *faces = &_cellIds[_nbDownCells * cellId];
2011 for (int i = 0; i < _nbDownCells; i++)
2015 faces[i] = lowCellId;
2018 if (faces[i] == lowCellId)
2022 // MESSAGE("-------------------------------------> trop de faces ! " << cellId << " " << lowCellId);
2025 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2026 * The hexahedron is defined by the eight points (0-7), where (0,1,2,3) is the base
2027 * of the hexahedron which, using the right hand rule, forms a quadrilateral whose normal
2028 * points in the direction of the opposite face (4,5,6,7).
2029 * @see vtkHexahedron.h in Filtering
2030 * @param cellId volumeId in vtkUnstructuredGrid
2031 * @param facesWithNodes vector of face descriptors to be filled
2033 void SMDS_DownHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2035 // --- find point id's of the volume
2038 vtkIdType *nodes; // will refer to the point id's of the volume
2039 _grid->GetCellPoints(cellId, npts, nodes);
2041 // --- create all the ordered list of node id's for each face
2043 facesWithNodes.nbElems = 6;
2045 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2046 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2047 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2048 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2049 facesWithNodes.elems[0].nbNodes = 4;
2050 facesWithNodes.elems[0].vtkType = VTK_QUAD;
2052 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2053 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2054 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2055 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2056 facesWithNodes.elems[1].nbNodes = 4;
2057 facesWithNodes.elems[1].vtkType = VTK_QUAD;
2059 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2060 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2061 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2062 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2063 facesWithNodes.elems[2].nbNodes = 4;
2064 facesWithNodes.elems[2].vtkType = VTK_QUAD;
2066 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2067 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2068 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2069 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2070 facesWithNodes.elems[3].nbNodes = 4;
2071 facesWithNodes.elems[3].vtkType = VTK_QUAD;
2073 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2074 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2075 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2076 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2077 facesWithNodes.elems[4].nbNodes = 4;
2078 facesWithNodes.elems[4].vtkType = VTK_QUAD;
2080 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2081 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2082 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2083 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2084 facesWithNodes.elems[5].nbNodes = 4;
2085 facesWithNodes.elems[5].vtkType = VTK_QUAD;
2088 // ---------------------------------------------------------------------------
2090 SMDS_DownQuadHexa::SMDS_DownQuadHexa(SMDS_UnstructuredGrid *grid) :
2091 SMDS_Down3D(grid, 6)
2093 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2094 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2095 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2096 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2097 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2098 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2101 SMDS_DownQuadHexa::~SMDS_DownQuadHexa()
2105 void SMDS_DownQuadHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
2109 for (int i = 0; i < orderedNodes.size(); i++)
2110 setNodes.insert(orderedNodes[i]);
2111 //MESSAGE("cellId = " << cellId);
2114 vtkIdType *nodes; // will refer to the point id's of the volume
2115 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
2118 //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};
2119 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,
2120 4, 0, 1, 5,16, 8,17,12, 5, 1, 2, 6,17, 9,18,13, 6, 2, 3, 7,18,10,19,14};
2121 for (int k = 0; k < 6; k++)
2124 for (int i = 0; i < 8; i++)
2125 tofind.insert(nodes[ids[8 * k + i]]);
2126 if (setNodes == tofind)
2128 for (int i = 0; i < 8; i++)
2129 orderedNodes[i] = nodes[ids[8 * k + i]];
2133 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
2134 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2135 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2138 void SMDS_DownQuadHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
2140 //ASSERT((cellId >=0)&& (cellId < _maxId));
2141 int *faces = &_cellIds[_nbDownCells * cellId];
2142 for (int i = 0; i < _nbDownCells; i++)
2146 faces[i] = lowCellId;
2149 if (faces[i] == lowCellId)
2155 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2156 * The ordering of the twenty points defining the quadratic hexahedron cell is point id's (0-7,8-19)
2157 * where point id's 0-7 are the eight corner vertices of the cube, followed by twelve mid-edge nodes (8-19).
2158 * Note that these mid-edge nodes lie on the edges defined by
2159 * (0,1), (1,2), (2,3), (3,0), (4,5), (5,6), (6,7), (7,4), (0,4), (1,5), (2,6), (3,7).
2160 * @see vtkQuadraticHexahedron.h in Filtering
2161 * @param cellId volumeId in vtkUnstructuredGrid
2162 * @param facesWithNodes vector of face descriptors to be filled
2164 void SMDS_DownQuadHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2166 // --- find point id's of the volume
2169 vtkIdType *nodes; // will refer to the point id's of the volume
2170 _grid->GetCellPoints(cellId, npts, nodes);
2172 // --- create all the ordered list of node id's for each face
2174 facesWithNodes.nbElems = 6;
2176 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2177 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2178 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2179 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2180 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
2181 facesWithNodes.elems[0].nodeIds[5] = nodes[9];
2182 facesWithNodes.elems[0].nodeIds[6] = nodes[10];
2183 facesWithNodes.elems[0].nodeIds[7] = nodes[11];
2184 facesWithNodes.elems[0].nbNodes = 8;
2185 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
2187 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2188 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2189 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2190 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2191 facesWithNodes.elems[1].nodeIds[4] = nodes[12];
2192 facesWithNodes.elems[1].nodeIds[5] = nodes[13];
2193 facesWithNodes.elems[1].nodeIds[6] = nodes[14];
2194 facesWithNodes.elems[1].nodeIds[7] = nodes[15];
2195 facesWithNodes.elems[1].nbNodes = 8;
2196 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
2198 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2199 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2200 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2201 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2202 facesWithNodes.elems[2].nodeIds[4] = nodes[8];
2203 facesWithNodes.elems[2].nodeIds[5] = nodes[17];
2204 facesWithNodes.elems[2].nodeIds[6] = nodes[12];
2205 facesWithNodes.elems[2].nodeIds[7] = nodes[16];
2206 facesWithNodes.elems[2].nbNodes = 8;
2207 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
2209 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2210 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2211 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2212 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2213 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
2214 facesWithNodes.elems[3].nodeIds[5] = nodes[18];
2215 facesWithNodes.elems[3].nodeIds[6] = nodes[13];
2216 facesWithNodes.elems[3].nodeIds[7] = nodes[17];
2217 facesWithNodes.elems[3].nbNodes = 8;
2218 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_QUAD;
2220 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2221 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2222 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2223 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2224 facesWithNodes.elems[4].nodeIds[4] = nodes[18];
2225 facesWithNodes.elems[4].nodeIds[5] = nodes[14];
2226 facesWithNodes.elems[4].nodeIds[6] = nodes[19];
2227 facesWithNodes.elems[4].nodeIds[7] = nodes[10];
2228 facesWithNodes.elems[4].nbNodes = 8;
2229 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_QUAD;
2231 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2232 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2233 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2234 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2235 facesWithNodes.elems[5].nodeIds[4] = nodes[19];
2236 facesWithNodes.elems[5].nodeIds[5] = nodes[15];
2237 facesWithNodes.elems[5].nodeIds[6] = nodes[16];
2238 facesWithNodes.elems[5].nodeIds[7] = nodes[11];
2239 facesWithNodes.elems[5].nbNodes = 8;
2240 facesWithNodes.elems[5].vtkType = VTK_QUADRATIC_QUAD;
2243 // ---------------------------------------------------------------------------