1 // Copyright (C) 2010-2011 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_TETRA] = 3;
56 _cellDimension[VTK_QUADRATIC_TETRA] = 3;
57 _cellDimension[VTK_HEXAHEDRON] = 3;
58 _cellDimension[VTK_QUADRATIC_HEXAHEDRON] = 3;
59 _cellDimension[VTK_WEDGE] = 3;
60 _cellDimension[VTK_QUADRATIC_WEDGE] = 3;
61 _cellDimension[VTK_PYRAMID] = 3;
62 _cellDimension[VTK_QUADRATIC_PYRAMID] = 3;
64 return _cellDimension[cellType];
67 // ---------------------------------------------------------------------------
69 /*! Generic constructor for all the downward connectivity structures (one per vtk cell type).
70 * The static structure for cell dimension is set only once.
71 * @param grid unstructured grid associated to the mesh.
72 * @param nbDownCells number of downward entities associated to this vtk type of cell.
75 SMDS_Downward::SMDS_Downward(SMDS_UnstructuredGrid *grid, int nbDownCells) :
76 _grid(grid), _nbDownCells(nbDownCells)
79 this->_cellIds.clear();
80 this->_cellTypes.clear();
81 if (_cellDimension.empty())
83 _cellDimension.resize(VTK_MAXTYPE + 1, 0);
84 _cellDimension[VTK_LINE] = 1;
85 _cellDimension[VTK_QUADRATIC_EDGE] = 1;
86 _cellDimension[VTK_TRIANGLE] = 2;
87 _cellDimension[VTK_QUADRATIC_TRIANGLE] = 2;
88 _cellDimension[VTK_QUAD] = 2;
89 _cellDimension[VTK_QUADRATIC_QUAD] = 2;
90 _cellDimension[VTK_TETRA] = 3;
91 _cellDimension[VTK_QUADRATIC_TETRA] = 3;
92 _cellDimension[VTK_HEXAHEDRON] = 3;
93 _cellDimension[VTK_QUADRATIC_HEXAHEDRON] = 3;
94 _cellDimension[VTK_WEDGE] = 3;
95 _cellDimension[VTK_QUADRATIC_WEDGE] = 3;
96 _cellDimension[VTK_PYRAMID] = 3;
97 _cellDimension[VTK_QUADRATIC_PYRAMID] = 3;
101 SMDS_Downward::~SMDS_Downward()
105 /*! Give or create an entry for downward connectivity structure relative to a cell.
106 * If the entry already exists, just return its id, otherwise, create it.
107 * The internal storage memory is allocated if needed.
108 * The SMDS_UnstructuredGrid::_cellIdToDownId vector is completed for vtkUnstructuredGrid cells.
109 * @param vtkId for a vtkUnstructuredGrid cell or -1 (default) for a created downward cell.
110 * @return the rank in downward[vtkType] structure.
112 int SMDS_Downward::addCell(int vtkId)
116 localId = _grid->CellIdToDownId(vtkId);
120 localId = this->_maxId;
122 this->allocate(_maxId);
125 this->_vtkCellIds[localId] = vtkId;
126 _grid->setCellIdToDownId(vtkId, localId);
128 this->initCell(localId);
132 /*! generic method do nothing. see derived methods
136 void SMDS_Downward::initCell(int cellId)
140 /*! Get the number of downward entities associated to a cell (always the same for a given vtk type of cell)
142 * @param cellId not used here.
145 int SMDS_Downward::getNumberOfDownCells(int cellId)
150 /*! get a pointer on the downward entities id's associated to a cell.
151 * @see SMDS_Downward::getNumberOfDownCells for the number of downward entities.
152 * @see SMDS_Downward::getDownTypes for the vtk cell types associated to the downward entities.
153 * @param cellId index of the cell in the downward structure relative to a given vtk cell type.
154 * @return table of downward entities id's.
156 const int* SMDS_Downward::getDownCells(int cellId)
158 //ASSERT((cellId >=0) && (cellId < _maxId));
159 return &_cellIds[_nbDownCells * cellId];
162 /*! get a list of vtk cell types associated to downward entities of a given cell, in the same order
163 * than the downward entities id's list (@see SMDS_Downward::getDownCells).
165 * @param cellId index of the cell in the downward structure relative to a vtk cell type.
166 * @return table of downward entities types.
168 const unsigned char* SMDS_Downward::getDownTypes(int cellId)
170 return &_cellTypes[0];
173 /*! add a downward entity of dimension n-1 (cell or node) to a given cell.
174 * Actual implementation is done in derived methods.
175 * @param cellId index of the parent cell (dimension n) in the downward structure relative to a vtk cell type.
176 * @param lowCellId index of the children cell to add (dimension n-1)
177 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
179 void SMDS_Downward::addDownCell(int cellId, int lowCellId, unsigned char aType)
181 ASSERT(0); // must be re-implemented in derived class
184 /*! add a downward entity of dimension n+1 to a given cell.
185 * Actual implementation is done in derived methods.
186 * @param cellId index of the children cell (dimension n) in the downward structure relative to a vtk cell type.
187 * @param upCellId index of the parent cell to add (dimension n+1)
188 * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
190 void SMDS_Downward::addUpCell(int cellId, int upCellId, unsigned char aType)
192 ASSERT(0); // must be re-implemented in derived class
195 int SMDS_Downward::getNodeSet(int cellId, int* nodeSet)
200 // ---------------------------------------------------------------------------
202 SMDS_Down1D::SMDS_Down1D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
203 SMDS_Downward(grid, nbDownCells)
205 _upCellIdsVector.clear();
206 _upCellTypesVector.clear();
208 _upCellTypes.clear();
209 _upCellIndex.clear();
212 SMDS_Down1D::~SMDS_Down1D()
216 /*! clear vectors used to reference 2D cells containing the edge
220 void SMDS_Down1D::initCell(int cellId)
222 _upCellIdsVector[cellId].clear();
223 _upCellTypesVector[cellId].clear();
226 /*! Resize the downward connectivity storage vector if needed.
228 * @param nbElems total number of elements of the same type required
230 void SMDS_Down1D::allocate(int nbElems)
232 if (nbElems >= _vtkCellIds.size())
234 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
235 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
236 _upCellIdsVector.resize(nbElems + SMDS_Mesh::chunkSize);
237 _upCellTypesVector.resize(nbElems + SMDS_Mesh::chunkSize);
241 void SMDS_Down1D::compactStorage()
243 _cellIds.resize(_nbDownCells * _maxId);
244 _vtkCellIds.resize(_maxId);
247 for (int i = 0; i < _maxId; i++)
248 sizeUpCells += _upCellIdsVector[i].size();
249 _upCellIds.resize(sizeUpCells, -1);
250 _upCellTypes.resize(sizeUpCells);
251 _upCellIndex.resize(_maxId + 1, -1); // id and types of rank i correspond to [ _upCellIndex[i], _upCellIndex[i+1] [
253 for (int i = 0; i < _maxId; i++)
255 _upCellIndex[i] = current;
256 for (int j = 0; j < _upCellIdsVector[i].size(); j++)
258 _upCellIds[current] = _upCellIdsVector[i][j];
259 _upCellTypes[current] = _upCellTypesVector[i][j];
263 _upCellIndex[_maxId] = current;
265 _upCellIdsVector.clear();
266 _upCellTypesVector.clear();
269 void SMDS_Down1D::addUpCell(int cellId, int upCellId, unsigned char aType)
271 //ASSERT((cellId >=0) && (cellId < _maxId));
272 int nbFaces = _upCellIdsVector[cellId].size();
273 for (int i = 0; i < nbFaces; i++)
275 if ((_upCellIdsVector[cellId][i] == upCellId) && (_upCellTypesVector[cellId][i] == aType))
277 return; // already done
280 _upCellIdsVector[cellId].push_back(upCellId);
281 _upCellTypesVector[cellId].push_back(aType);
284 int SMDS_Down1D::getNumberOfUpCells(int cellId)
286 //ASSERT((cellId >=0) && (cellId < _maxId));
287 return _upCellIndex[cellId + 1] - _upCellIndex[cellId];
290 const int* SMDS_Down1D::getUpCells(int cellId)
292 //ASSERT((cellId >=0) && (cellId < _maxId));
293 return &_upCellIds[_upCellIndex[cellId]];
296 const unsigned char* SMDS_Down1D::getUpTypes(int cellId)
298 //ASSERT((cellId >=0) && (cellId < _maxId));
299 return &_upCellTypes[_upCellIndex[cellId]];
302 void SMDS_Down1D::getNodeIds(int cellId, std::set<int>& nodeSet)
304 for (int i = 0; i < _nbDownCells; i++)
305 nodeSet.insert(_cellIds[_nbDownCells * cellId + i]);
308 int SMDS_Down1D::getNodeSet(int cellId, int* nodeSet)
310 for (int i = 0; i < _nbDownCells; i++)
311 nodeSet[i] = _cellIds[_nbDownCells * cellId + i];
315 void SMDS_Down1D::setNodes(int cellId, int vtkId)
318 vtkIdType *pts; // will refer to the point id's of the face
319 _grid->GetCellPoints(vtkId, npts, pts);
320 // MESSAGE(vtkId << " " << npts << " " << _nbDownCells);
321 //ASSERT(npts == _nbDownCells);
322 for (int i = 0; i < npts; i++)
324 _cellIds[_nbDownCells * cellId + i] = pts[i];
328 void SMDS_Down1D::setNodes(int cellId, const int* nodeIds)
330 //ASSERT(nodeIds.size() == _nbDownCells);
331 for (int i = 0; i < _nbDownCells; i++)
333 _cellIds[_nbDownCells * cellId + i] = nodeIds[i];
337 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
338 * We keep in the list the cells that contains all the nodes, we keep only volumes and faces.
339 * @param cellId id of the edge in the downward structure
340 * @param vtkIds vector of vtk id's
341 * @return number of vtk cells (size of vector)
343 int SMDS_Down1D::computeVtkCells(int cellId, std::vector<int>& vtkIds)
347 // --- find all the cells the points belong to, and how many of the points belong to a given cell
349 int *pts = &_cellIds[_nbDownCells * cellId];
350 int ncells = this->computeVtkCells(pts, vtkIds);
354 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
356 * @param pts list of points id's defining an edge
357 * @param vtkIds vector of vtk id's
358 * @return number of vtk cells (size of vector)
360 int SMDS_Down1D::computeVtkCells(int *pts, std::vector<int>& vtkIds)
363 // --- find all the cells the points belong to, and how many of the points belong to a given cell
368 for (int i = 0; i < _nbDownCells; i++)
370 vtkIdType point = pts[i];
371 int numCells = _grid->GetLinks()->GetNcells(point);
372 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
373 for (int j = 0; j < numCells; j++)
375 int vtkCellId = cells[j];
377 for (int k = 0; k < cnt; k++)
379 if (cellIds[k] == vtkCellId)
388 cellIds[cnt] = vtkCellId;
390 // TODO ASSERT(cnt<1000);
396 // --- find the face and volume cells: they contains all the points and are of type volume or face
399 for (int i = 0; i < cnt; i++)
401 if (cellCnt[i] == _nbDownCells)
403 int vtkElemId = cellIds[i];
404 int vtkType = _grid->GetCellType(vtkElemId);
405 if (SMDS_Downward::getCellDimension(vtkType) > 1)
407 vtkIds.push_back(vtkElemId);
416 /*! Build the list of downward faces from a list of vtk cells.
418 * @param cellId id of the edge in the downward structure
419 * @param vtkIds vector of vtk id's
420 * @param downFaces vector of face id's in downward structures
421 * @param downTypes vector of face types
422 * @return number of downward faces
424 int SMDS_Down1D::computeFaces(int cellId, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
426 int *pts = &_cellIds[_nbDownCells * cellId];
427 int nbFaces = this->computeFaces(pts, vtkIds, nbcells, downFaces, downTypes);
431 /*! Build the list of downward faces from a list of vtk cells.
433 * @param pts list of points id's defining an edge
434 * @param vtkIds vector of vtk id's
435 * @param downFaces vector of face id's in downward structures
436 * @param downTypes vector of face types
437 * @return number of downward faces
439 int SMDS_Down1D::computeFaces(int* pts, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
442 for (int i = 0; i < nbcells; i++)
444 int vtkId = vtkIds[i];
445 int vtkType = _grid->GetCellType(vtkId);
446 if (SMDS_Downward::getCellDimension(vtkType) == 2)
448 int faceId = _grid->CellIdToDownId(vtkId);
449 downFaces[cnt] = faceId;
450 downTypes[cnt] = vtkType;
453 else if (SMDS_Downward::getCellDimension(vtkType) == 3)
455 int volId = _grid->CellIdToDownId(vtkId);
456 SMDS_Downward * downvol = _grid->getDownArray(vtkType);
457 //const int *downIds = downvol->getDownCells(volId);
458 const unsigned char* downTypesVol = downvol->getDownTypes(volId);
459 int nbFaces = downvol->getNumberOfDownCells(volId);
460 const int* faceIds = downvol->getDownCells(volId);
461 for (int n = 0; n < nbFaces; n++)
463 SMDS_Down2D *downFace = static_cast<SMDS_Down2D*> (_grid->getDownArray(downTypesVol[n]));
464 bool isInFace = downFace->isInFace(faceIds[n], pts, _nbDownCells);
467 bool alreadySet = false;
468 for (int k = 0; k < cnt; k++)
469 if (faceIds[n] == downFaces[k])
476 downFaces[cnt] = faceIds[n];
477 downTypes[cnt] = downTypesVol[n];
487 // ---------------------------------------------------------------------------
489 SMDS_Down2D::SMDS_Down2D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
490 SMDS_Downward(grid, nbDownCells)
493 _upCellTypes.clear();
498 SMDS_Down2D::~SMDS_Down2D()
502 int SMDS_Down2D::getNumberOfUpCells(int cellId)
505 if (_upCellIds[2 * cellId] >= 0)
507 if (_upCellIds[2 * cellId + 1] >= 0)
512 const int* SMDS_Down2D::getUpCells(int cellId)
514 //ASSERT((cellId >=0) && (cellId < _maxId));
515 return &_upCellIds[2 * cellId];
518 const unsigned char* SMDS_Down2D::getUpTypes(int cellId)
520 //ASSERT((cellId >=0) && (cellId < _maxId));
521 return &_upCellTypes[2 * cellId];
524 void SMDS_Down2D::getNodeIds(int cellId, std::set<int>& nodeSet)
526 for (int i = 0; i < _nbDownCells; i++)
528 int downCellId = _cellIds[_nbDownCells * cellId + i];
529 unsigned char cellType = _cellTypes[i];
530 this->_grid->getDownArray(cellType)->getNodeIds(downCellId, nodeSet);
534 /*! Find in vtkUnstructuredGrid the volumes containing a face already stored in vtkUnstructuredGrid.
535 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
536 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
537 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
538 * @param cellId the face cell id in vkUnstructuredGrid
539 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
540 * @return number of volumes (0, 1 or 2)
542 int SMDS_Down2D::computeVolumeIds(int cellId, int* ids)
544 // --- find point id's of the face
547 vtkIdType *pts; // will refer to the point id's of the face
548 _grid->GetCellPoints(cellId, npts, pts);
550 for (int i = 0; i < npts; i++)
551 nodes.push_back(pts[i]);
552 int nvol = this->computeVolumeIdsFromNodesFace(&nodes[0], npts, ids);
556 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
557 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
558 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
559 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
561 * @param ids a couple of vtkId, initialized at -1 (no parent volume)
562 * @return number of volumes (0, 1 or 2)
564 int SMDS_Down2D::computeVolumeIds(ElemByNodesType& faceByNodes, int* ids)
566 int nvol = this->computeVolumeIdsFromNodesFace(&faceByNodes.nodeIds[0], faceByNodes.nbNodes, ids);
570 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
571 * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
572 * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
573 * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
574 * @param pts array of vtk node id's
575 * @param npts number of nodes
577 * @return number of volumes (0, 1 or 2)
579 int SMDS_Down2D::computeVolumeIdsFromNodesFace(int* pts, int npts, int* ids)
582 // --- find all the cells the points belong to, and how many of the points belong to a given cell
587 for (int i = 0; i < npts; i++)
589 vtkIdType point = pts[i];
590 int numCells = _grid->GetLinks()->GetNcells(point);
591 //MESSAGE("cells pour " << i << " " << numCells);
592 vtkIdType *cells = _grid->GetLinks()->GetCells(point);
593 for (int j = 0; j < numCells; j++)
595 int vtkCellId = cells[j];
597 for (int k = 0; k < cnt; k++)
599 if (cellIds[k] == vtkCellId)
608 cellIds[cnt] = vtkCellId;
610 // TODO ASSERT(cnt<1000);
616 // --- find the volume cells: they contains all the points and are of type volume
619 for (int i = 0; i < cnt; i++)
621 //MESSAGE("cell " << cellIds[i] << " points " << cellCnt[i]);
622 if (cellCnt[i] == npts)
624 int vtkElemId = cellIds[i];
625 int vtkType = _grid->GetCellType(vtkElemId);
626 if (SMDS_Downward::getCellDimension(vtkType) == 3)
628 ids[nvol] = vtkElemId; // store the volume id in given vector
639 void SMDS_Down2D::setTempNodes(int cellId, int vtkId)
642 vtkIdType *pts; // will refer to the point id's of the face
643 _grid->GetCellPoints(vtkId, npts, pts);
644 // MESSAGE(vtkId << " " << npts << " " << _nbNodes);
645 //ASSERT(npts == _nbNodes);
646 for (int i = 0; i < npts; i++)
648 _tempNodes[_nbNodes * cellId + i] = pts[i];
652 void SMDS_Down2D::setTempNodes(int cellId, ElemByNodesType& faceByNodes)
654 for (int i = 0; i < faceByNodes.nbNodes; i++)
655 _tempNodes[_nbNodes * cellId + i] = faceByNodes.nodeIds[i];
658 /*! Find if all the nodes belongs to the face.
660 * @param cellId the face cell Id
661 * @param nodeSet set of node id's to be found in the face list of nodes
664 bool SMDS_Down2D::isInFace(int cellId, int *pts, int npts)
667 int *nodes = &_tempNodes[_nbNodes * cellId];
668 for (int j = 0; j < npts; j++)
671 for (int i = 0; i < _nbNodes; i++)
673 if (nodes[i] == point)
680 return (nbFound == npts);
683 /*! Resize the downward connectivity storage vector if needed.
685 * @param nbElems total number of elements of the same type required
687 void SMDS_Down2D::allocate(int nbElems)
689 if (nbElems >= _vtkCellIds.size())
691 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
692 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
693 _upCellIds.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
694 _upCellTypes.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
695 _tempNodes.resize(_nbNodes * (nbElems + SMDS_Mesh::chunkSize), -1);
699 void SMDS_Down2D::compactStorage()
701 _cellIds.resize(_nbDownCells * _maxId);
702 _upCellIds.resize(2 * _maxId);
703 _upCellTypes.resize(2 * _maxId);
704 _vtkCellIds.resize(_maxId);
708 void SMDS_Down2D::addUpCell(int cellId, int upCellId, unsigned char aType)
710 //ASSERT((cellId >=0)&& (cellId < _maxId));
711 int *vols = &_upCellIds[2 * cellId];
712 unsigned char *types = &_upCellTypes[2 * cellId];
713 for (int i = 0; i < 2; i++)
717 vols[i] = upCellId; // use non affected volume
721 if ((vols[i] == upCellId) && (types[i] == aType)) // already done
727 int SMDS_Down2D::getNodeSet(int cellId, int* nodeSet)
729 for (int i = 0; i < _nbNodes; i++)
730 nodeSet[i] = _tempNodes[_nbNodes * cellId + i];
734 int SMDS_Down2D::FindEdgeByNodes(int cellId, ElemByNodesType& edgeByNodes)
736 int *edges = &_cellIds[_nbDownCells * cellId];
737 for (int i = 0; i < _nbDownCells; i++)
739 if ((edges[i] >= 0) && (edgeByNodes.vtkType == _cellTypes[i]))
742 int npts = this->_grid->getDownArray(edgeByNodes.vtkType)->getNodeSet(edges[i], nodeSet);
744 for (int j = 0; j < npts; j++)
746 int point = edgeByNodes.nodeIds[j];
748 for (int k = 0; k < npts; k++)
750 if (nodeSet[k] == point)
766 // ---------------------------------------------------------------------------
768 SMDS_Down3D::SMDS_Down3D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
769 SMDS_Downward(grid, nbDownCells)
773 SMDS_Down3D::~SMDS_Down3D()
777 void SMDS_Down3D::allocate(int nbElems)
779 if (nbElems >= _vtkCellIds.size())
781 _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
782 _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
786 void SMDS_Down3D::compactStorage()
788 // nothing to do, size was known before
791 int SMDS_Down3D::getNumberOfUpCells(int cellId)
796 const int* SMDS_Down3D::getUpCells(int cellId)
801 const unsigned char* SMDS_Down3D::getUpTypes(int cellId)
806 void SMDS_Down3D::getNodeIds(int cellId, std::set<int>& nodeSet)
808 int vtkId = this->_vtkCellIds[cellId];
810 vtkIdType *nodes; // will refer to the point id's of the volume
811 _grid->GetCellPoints(vtkId, npts, nodes);
812 for (int i = 0; i < npts; i++)
813 nodeSet.insert(nodes[i]);
816 int SMDS_Down3D::FindFaceByNodes(int cellId, ElemByNodesType& faceByNodes)
818 int *faces = &_cellIds[_nbDownCells * cellId];
821 for (int i = 0; i < _nbDownCells; i++)
823 if ((faces[i] >= 0) && (faceByNodes.vtkType == _cellTypes[i]))
826 npoints = faceByNodes.nbNodes;
829 int npts = this->_grid->getDownArray(faceByNodes.vtkType)->getNodeSet(faces[i], nodeSet);
831 continue; // skip this face
833 for (int j = 0; j < npts; j++)
835 int point = faceByNodes.nodeIds[j];
837 for (int k = 0; k < npts; k++)
839 if (nodeSet[k] == point)
842 break; // point j is in the 2 faces, skip remaining k values
846 break; // point j is not in the 2 faces, skip the remaining tests
855 // ---------------------------------------------------------------------------
857 SMDS_DownEdge::SMDS_DownEdge(SMDS_UnstructuredGrid *grid) :
860 _cellTypes.push_back(VTK_VERTEX);
861 _cellTypes.push_back(VTK_VERTEX);
864 SMDS_DownEdge::~SMDS_DownEdge()
868 // ---------------------------------------------------------------------------
870 SMDS_DownQuadEdge::SMDS_DownQuadEdge(SMDS_UnstructuredGrid *grid) :
873 _cellTypes.push_back(VTK_VERTEX);
874 _cellTypes.push_back(VTK_VERTEX);
875 _cellTypes.push_back(VTK_VERTEX);
878 SMDS_DownQuadEdge::~SMDS_DownQuadEdge()
882 // ---------------------------------------------------------------------------
884 SMDS_DownTriangle::SMDS_DownTriangle(SMDS_UnstructuredGrid *grid) :
887 _cellTypes.push_back(VTK_LINE);
888 _cellTypes.push_back(VTK_LINE);
889 _cellTypes.push_back(VTK_LINE);
893 SMDS_DownTriangle::~SMDS_DownTriangle()
897 void SMDS_DownTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
899 int *nodes = &_tempNodes[_nbNodes * cellId];
900 edgesWithNodes.nbElems = 3;
902 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
903 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
904 edgesWithNodes.elems[0].nbNodes = 2;
905 edgesWithNodes.elems[0].vtkType = VTK_LINE;
907 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
908 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
909 edgesWithNodes.elems[1].nbNodes = 2;
910 edgesWithNodes.elems[1].vtkType = VTK_LINE;
912 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
913 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
914 edgesWithNodes.elems[2].nbNodes = 2;
915 edgesWithNodes.elems[2].vtkType = VTK_LINE;
918 void SMDS_DownTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
920 //ASSERT((cellId >=0)&& (cellId < _maxId));
921 //ASSERT(aType == VTK_LINE);
922 int *faces = &_cellIds[_nbDownCells * cellId];
923 for (int i = 0; i < _nbDownCells; i++)
927 faces[i] = lowCellId;
930 if (faces[i] == lowCellId)
936 // ---------------------------------------------------------------------------
938 SMDS_DownQuadTriangle::SMDS_DownQuadTriangle(SMDS_UnstructuredGrid *grid) :
941 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
942 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
943 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
947 SMDS_DownQuadTriangle::~SMDS_DownQuadTriangle()
951 void SMDS_DownQuadTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
953 int *nodes = &_tempNodes[_nbNodes * cellId];
954 edgesWithNodes.nbElems = 3;
956 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
957 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
958 edgesWithNodes.elems[0].nodeIds[2] = nodes[3];
959 edgesWithNodes.elems[0].nbNodes = 3;
960 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
962 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
963 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
964 edgesWithNodes.elems[1].nodeIds[2] = nodes[4];
965 edgesWithNodes.elems[1].nbNodes = 3;
966 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
968 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
969 edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
970 edgesWithNodes.elems[2].nodeIds[2] = nodes[5];
971 edgesWithNodes.elems[2].nbNodes = 3;
972 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
975 void SMDS_DownQuadTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
977 //ASSERT((cellId >=0)&& (cellId < _maxId));
978 //ASSERT(aType == VTK_QUADRATIC_EDGE);
979 int *edges = &_cellIds[_nbDownCells * cellId];
980 for (int i = 0; i < _nbDownCells; i++)
984 edges[i] = lowCellId;
987 if (edges[i] == lowCellId)
993 // ---------------------------------------------------------------------------
995 SMDS_DownQuadrangle::SMDS_DownQuadrangle(SMDS_UnstructuredGrid *grid) :
998 _cellTypes.push_back(VTK_LINE);
999 _cellTypes.push_back(VTK_LINE);
1000 _cellTypes.push_back(VTK_LINE);
1001 _cellTypes.push_back(VTK_LINE);
1005 SMDS_DownQuadrangle::~SMDS_DownQuadrangle()
1009 void SMDS_DownQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
1011 int *nodes = &_tempNodes[_nbNodes * cellId];
1012 edgesWithNodes.nbElems = 4;
1014 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1015 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1016 edgesWithNodes.elems[0].nbNodes = 2;
1017 edgesWithNodes.elems[0].vtkType = VTK_LINE;
1019 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1020 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1021 edgesWithNodes.elems[1].nbNodes = 2;
1022 edgesWithNodes.elems[1].vtkType = VTK_LINE;
1024 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1025 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1026 edgesWithNodes.elems[2].nbNodes = 2;
1027 edgesWithNodes.elems[2].vtkType = VTK_LINE;
1029 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1030 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1031 edgesWithNodes.elems[3].nbNodes = 2;
1032 edgesWithNodes.elems[3].vtkType = VTK_LINE;
1035 void SMDS_DownQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1037 //ASSERT((cellId >=0)&& (cellId < _maxId));
1038 //ASSERT(aType == VTK_LINE);
1039 int *faces = &_cellIds[_nbDownCells * cellId];
1040 for (int i = 0; i < _nbDownCells; i++)
1044 faces[i] = lowCellId;
1047 if (faces[i] == lowCellId)
1053 // ---------------------------------------------------------------------------
1055 SMDS_DownQuadQuadrangle::SMDS_DownQuadQuadrangle(SMDS_UnstructuredGrid *grid) :
1056 SMDS_Down2D(grid, 4)
1058 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1059 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1060 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1061 _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1065 SMDS_DownQuadQuadrangle::~SMDS_DownQuadQuadrangle()
1069 void SMDS_DownQuadQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
1071 int *nodes = &_tempNodes[_nbNodes * cellId];
1072 edgesWithNodes.nbElems = 4;
1074 edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1075 edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1076 edgesWithNodes.elems[0].nodeIds[2] = nodes[4];
1077 edgesWithNodes.elems[0].nbNodes = 3;
1078 edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
1080 edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1081 edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1082 edgesWithNodes.elems[1].nodeIds[2] = nodes[5];
1083 edgesWithNodes.elems[1].nbNodes = 3;
1084 edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
1086 edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1087 edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1088 edgesWithNodes.elems[2].nodeIds[2] = nodes[6];
1089 edgesWithNodes.elems[2].nbNodes = 3;
1090 edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
1092 edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1093 edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1094 edgesWithNodes.elems[3].nodeIds[2] = nodes[7];
1095 edgesWithNodes.elems[3].nbNodes = 3;
1096 edgesWithNodes.elems[3].vtkType = VTK_QUADRATIC_EDGE;
1099 void SMDS_DownQuadQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1101 //ASSERT((cellId >=0)&& (cellId < _maxId));
1102 //ASSERT(aType == VTK_QUADRATIC_EDGE);
1103 int *faces = &_cellIds[_nbDownCells * cellId];
1104 for (int i = 0; i < _nbDownCells; i++)
1108 faces[i] = lowCellId;
1111 if (faces[i] == lowCellId)
1117 // ---------------------------------------------------------------------------
1119 SMDS_DownTetra::SMDS_DownTetra(SMDS_UnstructuredGrid *grid) :
1120 SMDS_Down3D(grid, 4)
1122 _cellTypes.push_back(VTK_TRIANGLE);
1123 _cellTypes.push_back(VTK_TRIANGLE);
1124 _cellTypes.push_back(VTK_TRIANGLE);
1125 _cellTypes.push_back(VTK_TRIANGLE);
1128 SMDS_DownTetra::~SMDS_DownTetra()
1132 void SMDS_DownTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1136 for (int i = 0; i < orderedNodes.size(); i++)
1137 setNodes.insert(orderedNodes[i]);
1138 //MESSAGE("cellId = " << cellId);
1141 vtkIdType *nodes; // will refer to the point id's of the volume
1142 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1145 int ids[12] = { 0, 1, 2, 0, 3, 1, 2, 3, 0, 1, 3, 2 };
1146 //int ids[12] = { 2, 1, 0, 1, 3, 0, 0, 3, 2, 2, 3, 1 };
1147 for (int k = 0; k < 4; k++)
1150 for (int i = 0; i < 3; i++)
1151 tofind.insert(nodes[ids[3 * k + i]]);
1152 if (setNodes == tofind)
1154 for (int i = 0; i < 3; i++)
1155 orderedNodes[i] = nodes[ids[3 * k + i]];
1159 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1160 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1161 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1164 void SMDS_DownTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1166 //ASSERT((cellId >=0)&& (cellId < _maxId));
1167 //ASSERT(aType == VTK_TRIANGLE);
1168 int *faces = &_cellIds[_nbDownCells * cellId];
1169 for (int i = 0; i < _nbDownCells; i++)
1173 faces[i] = lowCellId;
1176 if (faces[i] == lowCellId)
1182 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1183 * The linear tetrahedron is defined by four points.
1184 * @see vtkTetra.h in Filtering.
1185 * @param cellId volumeId in vtkUnstructuredGrid
1186 * @param facesWithNodes vector of face descriptors to be filled
1188 void SMDS_DownTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1190 // --- find point id's of the volume
1193 vtkIdType *nodes; // will refer to the point id's of the volume
1194 _grid->GetCellPoints(cellId, npts, nodes);
1196 // --- create all the ordered list of node id's for each face
1198 facesWithNodes.nbElems = 4;
1200 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1201 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1202 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1203 facesWithNodes.elems[0].nbNodes = 3;
1204 facesWithNodes.elems[0].vtkType = VTK_TRIANGLE;
1206 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1207 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1208 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1209 facesWithNodes.elems[1].nbNodes = 3;
1210 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1212 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1213 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1214 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1215 facesWithNodes.elems[2].nbNodes = 3;
1216 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1218 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1219 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1220 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1221 facesWithNodes.elems[3].nbNodes = 3;
1222 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1225 // ---------------------------------------------------------------------------
1227 SMDS_DownQuadTetra::SMDS_DownQuadTetra(SMDS_UnstructuredGrid *grid) :
1228 SMDS_Down3D(grid, 4)
1230 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1231 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1232 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1233 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1236 SMDS_DownQuadTetra::~SMDS_DownQuadTetra()
1240 void SMDS_DownQuadTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1244 for (int i = 0; i < orderedNodes.size(); i++)
1245 setNodes.insert(orderedNodes[i]);
1246 //MESSAGE("cellId = " << cellId);
1249 vtkIdType *nodes; // will refer to the point id's of the volume
1250 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1253 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 };
1254 //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 };
1255 for (int k = 0; k < 4; k++)
1258 for (int i = 0; i < 6; i++)
1259 tofind.insert(nodes[ids[6 * k + i]]);
1260 if (setNodes == tofind)
1262 for (int i = 0; i < 6; i++)
1263 orderedNodes[i] = nodes[ids[6 * k + i]];
1267 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1268 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1269 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1272 void SMDS_DownQuadTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1274 //ASSERT((cellId >=0)&& (cellId < _maxId));
1275 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1276 int *faces = &_cellIds[_nbDownCells * cellId];
1277 for (int i = 0; i < _nbDownCells; i++)
1281 faces[i] = lowCellId;
1284 if (faces[i] == lowCellId)
1290 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1291 * The ordering of the ten points defining the quadratic tetrahedron cell is point id's (0-3,4-9)
1292 * where id's 0-3 are the four tetrahedron vertices;
1293 * 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).
1294 * @see vtkQuadraticTetra.h in Filtering.
1295 * @param cellId volumeId in vtkUnstructuredGrid
1296 * @param facesWithNodes vector of face descriptors to be filled
1298 void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1300 // --- find point id's of the volume
1303 vtkIdType *nodes; // will refer to the point id's of the volume
1304 _grid->GetCellPoints(cellId, npts, nodes);
1306 // --- create all the ordered list of node id's for each face
1308 facesWithNodes.nbElems = 4;
1310 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1311 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1312 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1313 facesWithNodes.elems[0].nodeIds[3] = nodes[4];
1314 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1315 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1316 facesWithNodes.elems[0].nbNodes = 6;
1317 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_TRIANGLE;
1319 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1320 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1321 facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1322 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1323 facesWithNodes.elems[1].nodeIds[4] = nodes[8];
1324 facesWithNodes.elems[1].nodeIds[5] = nodes[7];
1325 facesWithNodes.elems[1].nbNodes = 6;
1326 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1328 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1329 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1330 facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1331 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1332 facesWithNodes.elems[2].nodeIds[4] = nodes[9];
1333 facesWithNodes.elems[2].nodeIds[5] = nodes[7];
1334 facesWithNodes.elems[2].nbNodes = 6;
1335 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1337 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1338 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1339 facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1340 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
1341 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
1342 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1343 facesWithNodes.elems[3].nbNodes = 6;
1344 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1347 // ---------------------------------------------------------------------------
1349 SMDS_DownPyramid::SMDS_DownPyramid(SMDS_UnstructuredGrid *grid) :
1350 SMDS_Down3D(grid, 5)
1352 _cellTypes.push_back(VTK_QUAD);
1353 _cellTypes.push_back(VTK_TRIANGLE);
1354 _cellTypes.push_back(VTK_TRIANGLE);
1355 _cellTypes.push_back(VTK_TRIANGLE);
1356 _cellTypes.push_back(VTK_TRIANGLE);
1359 SMDS_DownPyramid::~SMDS_DownPyramid()
1363 void SMDS_DownPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1367 for (int i = 0; i < orderedNodes.size(); i++)
1368 setNodes.insert(orderedNodes[i]);
1369 //MESSAGE("cellId = " << cellId);
1372 vtkIdType *nodes; // will refer to the point id's of the volume
1373 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1376 int ids[16] = { 0, 1, 2, 3, 0, 3, 4, 3, 2, 4, 2, 1, 4, 1, 0, 4 };
1379 for (int i = 0; i < 4; i++)
1380 tofind.insert(nodes[ids[i]]);
1381 if (setNodes == tofind)
1383 for (int i = 0; i < 4; i++)
1384 orderedNodes[i] = nodes[ids[i]];
1387 for (int k = 0; k < 4; k++)
1390 for (int i = 0; i < 3; i++)
1391 tofind.insert(nodes[ids[4 + 3 * k + i]]);
1392 if (setNodes == tofind)
1394 for (int i = 0; i < 3; i++)
1395 orderedNodes[i] = nodes[ids[4 + 3 * k + i]];
1399 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1400 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1401 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1404 void SMDS_DownPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1406 //ASSERT((cellId >=0) && (cellId < _maxId));
1407 int *faces = &_cellIds[_nbDownCells * cellId];
1408 if (aType == VTK_QUAD)
1412 faces[0] = lowCellId;
1415 if (faces[0] == lowCellId)
1420 //ASSERT(aType == VTK_TRIANGLE);
1421 for (int i = 1; i < _nbDownCells; i++)
1425 faces[i] = lowCellId;
1428 if (faces[i] == lowCellId)
1435 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1436 * The pyramid is defined by the five points (0-4) where (0,1,2,3) is the base of the pyramid which,
1437 * using the right hand rule, forms a quadrilateral whose normal points in the direction of the
1438 * pyramid apex at vertex #4.
1439 * @see vtkPyramid.h in Filtering.
1440 * @param cellId volumeId in vtkUnstructuredGrid
1441 * @param facesWithNodes vector of face descriptors to be filled
1443 void SMDS_DownPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1445 // --- find point id's of the volume
1448 vtkIdType *nodes; // will refer to the point id's of the volume
1449 _grid->GetCellPoints(cellId, npts, nodes);
1451 // --- create all the ordered list of node id's for each face
1453 facesWithNodes.nbElems = 5;
1455 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1456 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1457 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1458 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1459 facesWithNodes.elems[0].nbNodes = 4;
1460 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1462 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1463 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1464 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1465 facesWithNodes.elems[1].nbNodes = 3;
1466 facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1468 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1469 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1470 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1471 facesWithNodes.elems[2].nbNodes = 3;
1472 facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1474 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1475 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1476 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1477 facesWithNodes.elems[3].nbNodes = 3;
1478 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1480 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1481 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1482 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1483 facesWithNodes.elems[4].nbNodes = 3;
1484 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1487 // ---------------------------------------------------------------------------
1489 SMDS_DownQuadPyramid::SMDS_DownQuadPyramid(SMDS_UnstructuredGrid *grid) :
1490 SMDS_Down3D(grid, 5)
1492 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1493 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1494 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1495 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1496 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1499 SMDS_DownQuadPyramid::~SMDS_DownQuadPyramid()
1503 void SMDS_DownQuadPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1507 for (int i = 0; i < orderedNodes.size(); i++)
1508 setNodes.insert(orderedNodes[i]);
1509 //MESSAGE("cellId = " << cellId);
1512 vtkIdType *nodes; // will refer to the point id's of the volume
1513 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1516 int ids[32] = { 0, 1, 2, 3, 5, 6, 7, 8,
1517 0, 3, 4, 8, 12, 9, 3, 2, 4, 7 , 11, 12, 2, 1, 4, 6, 10, 11, 1, 0, 4, 5, 9, 10 };
1520 for (int i = 0; i < 4; i++)
1521 tofind.insert(nodes[ids[i]]);
1522 if (setNodes == tofind)
1524 for (int i = 0; i < 8; i++)
1525 orderedNodes[i] = nodes[ids[i]];
1528 for (int k = 0; k < 4; k++)
1531 for (int i = 0; i < 6; i++)
1532 tofind.insert(nodes[ids[8 + 6 * k + i]]);
1533 if (setNodes == tofind)
1535 for (int i = 0; i < 6; i++)
1536 orderedNodes[i] = nodes[ids[8 + 6 * k + i]];
1540 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1541 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1542 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1545 void SMDS_DownQuadPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1547 //ASSERT((cellId >=0) && (cellId < _maxId));
1548 int *faces = &_cellIds[_nbDownCells * cellId];
1549 if (aType == VTK_QUADRATIC_QUAD)
1553 faces[0] = lowCellId;
1556 if (faces[0] == lowCellId)
1561 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1562 for (int i = 1; i < _nbDownCells; i++)
1566 faces[i] = lowCellId;
1569 if (faces[i] == lowCellId)
1576 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1577 * The ordering of the thirteen points defining the quadratic pyramid cell is point id's (0-4,5-12)
1578 * where point id's 0-4 are the five corner vertices of the pyramid; followed
1579 * by eight mid-edge nodes (5-12). Note that these mid-edge nodes lie on the edges defined by
1580 * 5(0,1), 6(1,2), 7(2,3), 8(3,0), 9(0,4), 10(1,4), 11(2,4), 12(3,4).
1581 * @see vtkQuadraticPyramid.h in Filtering.
1582 * @param cellId volumeId in vtkUnstructuredGrid
1583 * @param facesWithNodes vector of face descriptors to be filled
1585 void SMDS_DownQuadPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1587 // --- find point id's of the volume
1590 vtkIdType *nodes; // will refer to the point id's of the volume
1591 _grid->GetCellPoints(cellId, npts, nodes);
1593 // --- create all the ordered list of node id's for each face
1595 facesWithNodes.nbElems = 5;
1597 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1598 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1599 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1600 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1601 facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1602 facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1603 facesWithNodes.elems[0].nodeIds[6] = nodes[7];
1604 facesWithNodes.elems[0].nodeIds[7] = nodes[8];
1605 facesWithNodes.elems[0].nbNodes = 8;
1606 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1608 facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1609 facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1610 facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1611 facesWithNodes.elems[1].nodeIds[3] = nodes[5];
1612 facesWithNodes.elems[1].nodeIds[4] = nodes[10];
1613 facesWithNodes.elems[1].nodeIds[5] = nodes[9];
1614 facesWithNodes.elems[1].nbNodes = 6;
1615 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1617 facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1618 facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1619 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1620 facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1621 facesWithNodes.elems[2].nodeIds[4] = nodes[11];
1622 facesWithNodes.elems[2].nodeIds[5] = nodes[10];
1623 facesWithNodes.elems[2].nbNodes = 6;
1624 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1626 facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1627 facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1628 facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1629 facesWithNodes.elems[3].nodeIds[3] = nodes[7];
1630 facesWithNodes.elems[3].nodeIds[4] = nodes[12];
1631 facesWithNodes.elems[3].nodeIds[5] = nodes[11];
1632 facesWithNodes.elems[3].nbNodes = 6;
1633 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1635 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1636 facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1637 facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1638 facesWithNodes.elems[4].nodeIds[3] = nodes[8];
1639 facesWithNodes.elems[4].nodeIds[4] = nodes[9];
1640 facesWithNodes.elems[4].nodeIds[5] = nodes[12];
1641 facesWithNodes.elems[4].nbNodes = 6;
1642 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1645 // ---------------------------------------------------------------------------
1647 SMDS_DownPenta::SMDS_DownPenta(SMDS_UnstructuredGrid *grid) :
1648 SMDS_Down3D(grid, 5)
1650 _cellTypes.push_back(VTK_QUAD);
1651 _cellTypes.push_back(VTK_QUAD);
1652 _cellTypes.push_back(VTK_QUAD);
1653 _cellTypes.push_back(VTK_TRIANGLE);
1654 _cellTypes.push_back(VTK_TRIANGLE);
1657 SMDS_DownPenta::~SMDS_DownPenta()
1661 void SMDS_DownPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1665 for (int i = 0; i < orderedNodes.size(); i++)
1666 setNodes.insert(orderedNodes[i]);
1667 //MESSAGE("cellId = " << cellId);
1670 vtkIdType *nodes; // will refer to the point id's of the volume
1671 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1674 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1675 int ids[18] = { 0, 1, 2, 3, 5, 4, 0, 3, 4, 1, 1, 4, 5, 2, 2, 5, 3, 0 };
1677 for (int k = 0; k < 2; k++)
1680 for (int i = 0; i < 3; i++)
1681 tofind.insert(nodes[ids[3 * k + i]]);
1682 if (setNodes == tofind)
1684 for (int i = 0; i < 3; i++)
1685 orderedNodes[i] = nodes[ids[3 * k + i]];
1689 for (int k = 0; k < 3; k++)
1692 for (int i = 0; i < 4; i++)
1693 tofind.insert(nodes[ids[6 + 4 * k + i]]);
1694 if (setNodes == tofind)
1696 for (int i = 0; i < 4; i++)
1697 orderedNodes[i] = nodes[ids[6 + 4 * k + i]];
1701 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1702 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1703 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1706 void SMDS_DownPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1708 //ASSERT((cellId >=0) && (cellId < _maxId));
1709 int *faces = &_cellIds[_nbDownCells * cellId];
1710 if (aType == VTK_QUAD)
1711 for (int i = 0; i < 3; i++)
1715 faces[i] = lowCellId;
1718 if (faces[i] == lowCellId)
1723 //ASSERT(aType == VTK_TRIANGLE);
1724 for (int i = 3; i < _nbDownCells; i++)
1728 faces[i] = lowCellId;
1731 if (faces[i] == lowCellId)
1738 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's.
1739 * A wedge or pentahedron consists of two triangular and three quadrilateral faces
1740 * and is defined by the six points (0-5) where (0,1,2) is the base of the wedge which,
1741 * using the right hand rule, forms a triangle whose normal points outward
1742 * (away from the triangular face (3,4,5)).
1743 * @see vtkWedge.h in Filtering
1744 * @param cellId volumeId in vtkUnstructuredGrid
1745 * @param facesWithNodes vector of face descriptors to be filled
1747 void SMDS_DownPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1749 // --- find point id's of the volume
1752 vtkIdType *nodes; // will refer to the point id's of the volume
1753 _grid->GetCellPoints(cellId, npts, nodes);
1755 // --- create all the ordered list of node id's for each face
1757 facesWithNodes.nbElems = 5;
1759 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1760 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1761 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1762 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1763 facesWithNodes.elems[0].nbNodes = 4;
1764 facesWithNodes.elems[0].vtkType = VTK_QUAD;
1766 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1767 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1768 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1769 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1770 facesWithNodes.elems[1].nbNodes = 4;
1771 facesWithNodes.elems[1].vtkType = VTK_QUAD;
1773 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1774 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1775 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1776 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1777 facesWithNodes.elems[2].nbNodes = 4;
1778 facesWithNodes.elems[2].vtkType = VTK_QUAD;
1780 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1781 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1782 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1783 facesWithNodes.elems[3].nbNodes = 3;
1784 facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1786 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1787 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1788 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1789 facesWithNodes.elems[4].nbNodes = 3;
1790 facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1793 // ---------------------------------------------------------------------------
1795 SMDS_DownQuadPenta::SMDS_DownQuadPenta(SMDS_UnstructuredGrid *grid) :
1796 SMDS_Down3D(grid, 5)
1798 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1799 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1800 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1801 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1802 _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1805 SMDS_DownQuadPenta::~SMDS_DownQuadPenta()
1809 void SMDS_DownQuadPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1813 for (int i = 0; i < orderedNodes.size(); i++)
1814 setNodes.insert(orderedNodes[i]);
1815 //MESSAGE("cellId = " << cellId);
1818 vtkIdType *nodes; // will refer to the point id's of the volume
1819 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1822 //int ids[18] = { 0, 2, 1, 3, 4, 5, 0, 1, 4, 3, 1, 2, 5, 4, 2, 0, 3, 5 };
1823 int ids[36] = { 0, 1, 2, 6, 7, 8, 3, 5, 4, 11, 10, 9,
1824 0, 3, 4, 1, 12, 9, 13, 6, 1, 4, 5, 2, 13, 10, 14, 7, 2, 5, 3, 0, 14, 11, 12, 8 };
1826 for (int k = 0; k < 2; k++)
1829 for (int i = 0; i < 6; i++)
1830 tofind.insert(nodes[ids[6 * k + i]]);
1831 if (setNodes == tofind)
1833 for (int i = 0; i < 6; i++)
1834 orderedNodes[i] = nodes[ids[6 * k + i]];
1838 for (int k = 0; k < 3; k++)
1841 for (int i = 0; i < 8; i++)
1842 tofind.insert(nodes[ids[12 + 8 * k + i]]);
1843 if (setNodes == tofind)
1845 for (int i = 0; i < 8; i++)
1846 orderedNodes[i] = nodes[ids[12 + 8 * k + i]];
1850 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1851 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1852 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1855 void SMDS_DownQuadPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1857 //ASSERT((cellId >=0) && (cellId < _maxId));
1858 int *faces = &_cellIds[_nbDownCells * cellId];
1859 if (aType == VTK_QUADRATIC_QUAD)
1860 for (int i = 0; i < 3; i++)
1864 faces[i] = lowCellId;
1867 if (faces[i] == lowCellId)
1872 //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1873 for (int i = 3; i < _nbDownCells; i++)
1877 faces[i] = lowCellId;
1880 if (faces[i] == lowCellId)
1887 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
1888 * The quadratic wedge (or pentahedron) is defined by fifteen points.
1889 * The ordering of the fifteen points defining the cell is point id's (0-5,6-14)
1890 * where point id's 0-5 are the six corner vertices of the wedge, followed by
1891 * nine mid-edge nodes (6-14). Note that these mid-edge nodes lie on the edges defined by
1892 * (0,1), (1,2), (2,0), (3,4), (4,5), (5,3), (0,3), (1,4), (2,5).
1893 * @see vtkQuadraticWedge.h in Filtering
1894 * @param cellId volumeId in vtkUnstructuredGrid
1895 * @param facesWithNodes vector of face descriptors to be filled
1897 void SMDS_DownQuadPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1899 // --- find point id's of the volume
1902 vtkIdType *nodes; // will refer to the point id's of the volume
1903 _grid->GetCellPoints(cellId, npts, nodes);
1905 // --- create all the ordered list of node id's for each face
1907 facesWithNodes.nbElems = 5;
1909 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1910 facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1911 facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1912 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1913 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
1914 facesWithNodes.elems[0].nodeIds[5] = nodes[14];
1915 facesWithNodes.elems[0].nodeIds[6] = nodes[11];
1916 facesWithNodes.elems[0].nodeIds[7] = nodes[12];
1917 facesWithNodes.elems[0].nbNodes = 8;
1918 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1920 facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1921 facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1922 facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1923 facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1924 facesWithNodes.elems[1].nodeIds[4] = nodes[7];
1925 facesWithNodes.elems[1].nodeIds[5] = nodes[14];
1926 facesWithNodes.elems[1].nodeIds[6] = nodes[10];
1927 facesWithNodes.elems[1].nodeIds[7] = nodes[13];
1928 facesWithNodes.elems[1].nbNodes = 8;
1929 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
1931 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1932 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1933 facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1934 facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1935 facesWithNodes.elems[2].nodeIds[4] = nodes[6];
1936 facesWithNodes.elems[2].nodeIds[5] = nodes[13];
1937 facesWithNodes.elems[2].nodeIds[6] = nodes[9];
1938 facesWithNodes.elems[2].nodeIds[7] = nodes[12];
1939 facesWithNodes.elems[2].nbNodes = 8;
1940 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
1942 facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1943 facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1944 facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1945 facesWithNodes.elems[3].nodeIds[3] = nodes[6];
1946 facesWithNodes.elems[3].nodeIds[4] = nodes[7];
1947 facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1948 facesWithNodes.elems[3].nbNodes = 6;
1949 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1951 facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1952 facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1953 facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1954 facesWithNodes.elems[4].nodeIds[3] = nodes[9];
1955 facesWithNodes.elems[4].nodeIds[4] = nodes[10];
1956 facesWithNodes.elems[4].nodeIds[5] = nodes[11];
1957 facesWithNodes.elems[4].nbNodes = 6;
1958 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1961 // ---------------------------------------------------------------------------
1963 SMDS_DownHexa::SMDS_DownHexa(SMDS_UnstructuredGrid *grid) :
1964 SMDS_Down3D(grid, 6)
1966 _cellTypes.push_back(VTK_QUAD);
1967 _cellTypes.push_back(VTK_QUAD);
1968 _cellTypes.push_back(VTK_QUAD);
1969 _cellTypes.push_back(VTK_QUAD);
1970 _cellTypes.push_back(VTK_QUAD);
1971 _cellTypes.push_back(VTK_QUAD);
1974 SMDS_DownHexa::~SMDS_DownHexa()
1978 void SMDS_DownHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1982 for (int i = 0; i < orderedNodes.size(); i++)
1983 setNodes.insert(orderedNodes[i]);
1984 //MESSAGE("cellId = " << cellId);
1987 vtkIdType *nodes; // will refer to the point id's of the volume
1988 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1991 //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};
1992 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};
1993 for (int k = 0; k < 6; k++) // loop on the 6 faces
1996 for (int i = 0; i < 4; i++)
1997 tofind.insert(nodes[ids[4 * k + i]]); // node ids of the face i
1998 if (setNodes == tofind)
2000 for (int i = 0; i < 4; i++)
2001 orderedNodes[i] = nodes[ids[4 * k + i]];
2005 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
2006 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2007 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2008 MESSAGE(nodes[4] << " " << nodes[5] << " " << nodes[6] << " " << nodes[7]);
2011 void SMDS_DownHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
2013 //ASSERT((cellId >=0)&& (cellId < _maxId));
2014 int *faces = &_cellIds[_nbDownCells * cellId];
2015 for (int i = 0; i < _nbDownCells; i++)
2019 faces[i] = lowCellId;
2022 if (faces[i] == lowCellId)
2026 // MESSAGE("-------------------------------------> trop de faces ! " << cellId << " " << lowCellId);
2029 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2030 * The hexahedron is defined by the eight points (0-7), where (0,1,2,3) is the base
2031 * of the hexahedron which, using the right hand rule, forms a quadrilateral whose normal
2032 * points in the direction of the opposite face (4,5,6,7).
2033 * @see vtkHexahedron.h in Filtering
2034 * @param cellId volumeId in vtkUnstructuredGrid
2035 * @param facesWithNodes vector of face descriptors to be filled
2037 void SMDS_DownHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2039 // --- find point id's of the volume
2042 vtkIdType *nodes; // will refer to the point id's of the volume
2043 _grid->GetCellPoints(cellId, npts, nodes);
2045 // --- create all the ordered list of node id's for each face
2047 facesWithNodes.nbElems = 6;
2049 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2050 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2051 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2052 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2053 facesWithNodes.elems[0].nbNodes = 4;
2054 facesWithNodes.elems[0].vtkType = VTK_QUAD;
2056 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2057 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2058 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2059 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2060 facesWithNodes.elems[1].nbNodes = 4;
2061 facesWithNodes.elems[1].vtkType = VTK_QUAD;
2063 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2064 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2065 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2066 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2067 facesWithNodes.elems[2].nbNodes = 4;
2068 facesWithNodes.elems[2].vtkType = VTK_QUAD;
2070 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2071 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2072 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2073 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2074 facesWithNodes.elems[3].nbNodes = 4;
2075 facesWithNodes.elems[3].vtkType = VTK_QUAD;
2077 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2078 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2079 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2080 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2081 facesWithNodes.elems[4].nbNodes = 4;
2082 facesWithNodes.elems[4].vtkType = VTK_QUAD;
2084 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2085 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2086 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2087 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2088 facesWithNodes.elems[5].nbNodes = 4;
2089 facesWithNodes.elems[5].vtkType = VTK_QUAD;
2092 // ---------------------------------------------------------------------------
2094 SMDS_DownQuadHexa::SMDS_DownQuadHexa(SMDS_UnstructuredGrid *grid) :
2095 SMDS_Down3D(grid, 6)
2097 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2098 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2099 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2100 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2101 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2102 _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2105 SMDS_DownQuadHexa::~SMDS_DownQuadHexa()
2109 void SMDS_DownQuadHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
2113 for (int i = 0; i < orderedNodes.size(); i++)
2114 setNodes.insert(orderedNodes[i]);
2115 //MESSAGE("cellId = " << cellId);
2118 vtkIdType *nodes; // will refer to the point id's of the volume
2119 _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
2122 //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};
2123 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,
2124 4, 0, 1, 5,16, 8,17,12, 5, 1, 2, 6,17, 9,18,13, 6, 2, 3, 7,18,10,19,14};
2125 for (int k = 0; k < 6; k++)
2128 for (int i = 0; i < 8; i++)
2129 tofind.insert(nodes[ids[8 * k + i]]);
2130 if (setNodes == tofind)
2132 for (int i = 0; i < 8; i++)
2133 orderedNodes[i] = nodes[ids[8 * k + i]];
2137 MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
2138 MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2139 MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2142 void SMDS_DownQuadHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
2144 //ASSERT((cellId >=0)&& (cellId < _maxId));
2145 int *faces = &_cellIds[_nbDownCells * cellId];
2146 for (int i = 0; i < _nbDownCells; i++)
2150 faces[i] = lowCellId;
2153 if (faces[i] == lowCellId)
2159 /*! Create a list of faces described by a vtk Type and an ordered set of Node Id's
2160 * The ordering of the twenty points defining the quadratic hexahedron cell is point id's (0-7,8-19)
2161 * where point id's 0-7 are the eight corner vertices of the cube, followed by twelve mid-edge nodes (8-19).
2162 * Note that these mid-edge nodes lie on the edges defined by
2163 * (0,1), (1,2), (2,3), (3,0), (4,5), (5,6), (6,7), (7,4), (0,4), (1,5), (2,6), (3,7).
2164 * @see vtkQuadraticHexahedron.h in Filtering
2165 * @param cellId volumeId in vtkUnstructuredGrid
2166 * @param facesWithNodes vector of face descriptors to be filled
2168 void SMDS_DownQuadHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2170 // --- find point id's of the volume
2173 vtkIdType *nodes; // will refer to the point id's of the volume
2174 _grid->GetCellPoints(cellId, npts, nodes);
2176 // --- create all the ordered list of node id's for each face
2178 facesWithNodes.nbElems = 6;
2180 facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2181 facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2182 facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2183 facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2184 facesWithNodes.elems[0].nodeIds[4] = nodes[8];
2185 facesWithNodes.elems[0].nodeIds[5] = nodes[9];
2186 facesWithNodes.elems[0].nodeIds[6] = nodes[10];
2187 facesWithNodes.elems[0].nodeIds[7] = nodes[11];
2188 facesWithNodes.elems[0].nbNodes = 8;
2189 facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
2191 facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2192 facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2193 facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2194 facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2195 facesWithNodes.elems[1].nodeIds[4] = nodes[12];
2196 facesWithNodes.elems[1].nodeIds[5] = nodes[13];
2197 facesWithNodes.elems[1].nodeIds[6] = nodes[14];
2198 facesWithNodes.elems[1].nodeIds[7] = nodes[15];
2199 facesWithNodes.elems[1].nbNodes = 8;
2200 facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
2202 facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2203 facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2204 facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2205 facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2206 facesWithNodes.elems[2].nodeIds[4] = nodes[8];
2207 facesWithNodes.elems[2].nodeIds[5] = nodes[17];
2208 facesWithNodes.elems[2].nodeIds[6] = nodes[12];
2209 facesWithNodes.elems[2].nodeIds[7] = nodes[16];
2210 facesWithNodes.elems[2].nbNodes = 8;
2211 facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
2213 facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2214 facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2215 facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2216 facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2217 facesWithNodes.elems[3].nodeIds[4] = nodes[9];
2218 facesWithNodes.elems[3].nodeIds[5] = nodes[18];
2219 facesWithNodes.elems[3].nodeIds[6] = nodes[13];
2220 facesWithNodes.elems[3].nodeIds[7] = nodes[17];
2221 facesWithNodes.elems[3].nbNodes = 8;
2222 facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_QUAD;
2224 facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2225 facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2226 facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2227 facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2228 facesWithNodes.elems[4].nodeIds[4] = nodes[18];
2229 facesWithNodes.elems[4].nodeIds[5] = nodes[14];
2230 facesWithNodes.elems[4].nodeIds[6] = nodes[19];
2231 facesWithNodes.elems[4].nodeIds[7] = nodes[10];
2232 facesWithNodes.elems[4].nbNodes = 8;
2233 facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_QUAD;
2235 facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2236 facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2237 facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2238 facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2239 facesWithNodes.elems[5].nodeIds[4] = nodes[19];
2240 facesWithNodes.elems[5].nodeIds[5] = nodes[15];
2241 facesWithNodes.elems[5].nodeIds[6] = nodes[16];
2242 facesWithNodes.elems[5].nodeIds[7] = nodes[11];
2243 facesWithNodes.elems[5].nbNodes = 8;
2244 facesWithNodes.elems[5].vtkType = VTK_QUADRATIC_QUAD;
2247 // ---------------------------------------------------------------------------