]> SALOME platform Git repositories - modules/smesh.git/blob - src/SMDS/SMDS_Downward.cxx
Salome HOME
Correction of mistakes: 0022537: EDF 2877 SMESH : wrong information about the numbre...
[modules/smesh.git] / src / SMDS / SMDS_Downward.cxx
1 // Copyright (C) 2010-2014  CEA/DEN, EDF R&D, OPEN CASCADE
2 //
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
7 //
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.
12 //
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
16 //
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
18 //
19
20 // File: SMDS_Downward.cxx
21 // Created: Jun 3, 2010
22 // Author: prascle
23
24 #include "SMDS_Downward.hxx"
25 #include "SMDS_Mesh.hxx"
26 #include "utilities.h"
27
28 #include <vtkCellType.h>
29 #include <vtkCellLinks.h>
30
31 #include <map>
32
33 using namespace std;
34
35 // ---------------------------------------------------------------------------
36
37 vector<int> SMDS_Downward::_cellDimension;
38
39 /*! get the dimension of a cell (1,2,3 for 1D, 2D 3D) given the vtk cell type
40  *
41  * @param cellType vtk cell type @see vtkCellType.h
42  * @return 1,2 or 3
43  */
44 int SMDS_Downward::getCellDimension(unsigned char cellType)
45 {
46   if (_cellDimension.empty())
47     {
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_BIQUADRATIC_TRIANGLE] = 2;
54       _cellDimension[VTK_QUAD] = 2;
55       _cellDimension[VTK_QUADRATIC_QUAD] = 2;
56       _cellDimension[VTK_BIQUADRATIC_QUAD] = 2;
57       _cellDimension[VTK_TETRA] = 3;
58       _cellDimension[VTK_QUADRATIC_TETRA] = 3;
59       _cellDimension[VTK_HEXAHEDRON] = 3;
60       _cellDimension[VTK_QUADRATIC_HEXAHEDRON] = 3;
61       _cellDimension[VTK_TRIQUADRATIC_HEXAHEDRON] = 3;
62       _cellDimension[VTK_WEDGE] = 3;
63       _cellDimension[VTK_QUADRATIC_WEDGE] = 3;
64       _cellDimension[VTK_PYRAMID] = 3;
65       _cellDimension[VTK_QUADRATIC_PYRAMID] = 3;
66       _cellDimension[VTK_HEXAGONAL_PRISM] = 3;
67     }
68   return _cellDimension[cellType];
69 }
70
71 // ---------------------------------------------------------------------------
72
73 /*! Generic constructor for all the downward connectivity structures (one per vtk cell type).
74  *  The static structure for cell dimension is set only once.
75  * @param grid unstructured grid associated to the mesh.
76  * @param nbDownCells number of downward entities associated to this vtk type of cell.
77  * @return
78  */
79 SMDS_Downward::SMDS_Downward(SMDS_UnstructuredGrid *grid, int nbDownCells) :
80   _grid(grid), _nbDownCells(nbDownCells)
81 {
82   this->_maxId = 0;
83   this->_cellIds.clear();
84   this->_cellTypes.clear();
85   if (_cellDimension.empty())
86     getCellDimension( VTK_LINE );
87 }
88
89 SMDS_Downward::~SMDS_Downward()
90 {
91 }
92
93 /*! Give or create an entry for downward connectivity structure relative to a cell.
94  * If the entry already exists, just return its id, otherwise, create it.
95  * The internal storage memory is allocated if needed.
96  * The SMDS_UnstructuredGrid::_cellIdToDownId vector is completed for vtkUnstructuredGrid cells.
97  * @param vtkId for a vtkUnstructuredGrid cell  or -1 (default) for a created downward cell.
98  * @return the rank in downward[vtkType] structure.
99  */
100 int SMDS_Downward::addCell(int vtkId)
101 {
102   int localId = -1;
103   if (vtkId >= 0)
104     localId = _grid->CellIdToDownId(vtkId);
105   if (localId >= 0)
106     return localId;
107
108   localId = this->_maxId;
109   this->_maxId++;
110   this->allocate(_maxId);
111   if (vtkId >= 0)
112     {
113       this->_vtkCellIds[localId] = vtkId;
114       _grid->setCellIdToDownId(vtkId, localId);
115     }
116   this->initCell(localId);
117   return localId;
118 }
119
120 /*! generic method do nothing. see derived methods
121  *
122  * @param cellId
123  */
124 void SMDS_Downward::initCell(int cellId)
125 {
126 }
127
128 /*! Get the number of downward entities associated to a cell (always the same for a given vtk type of cell)
129  *
130  * @param cellId not used here.
131  * @return
132  */
133 int SMDS_Downward::getNumberOfDownCells(int cellId)
134 {
135   return _nbDownCells;
136 }
137
138 /*! get a pointer on the downward entities id's associated to a cell.
139  * @see SMDS_Downward::getNumberOfDownCells for the number of downward entities.
140  * @see SMDS_Downward::getDownTypes for the vtk cell types associated to the downward entities.
141  * @param cellId index of the cell in the downward structure relative to a given vtk cell type.
142  * @return table of downward entities id's.
143  */
144 const int* SMDS_Downward::getDownCells(int cellId)
145 {
146   //ASSERT((cellId >=0) && (cellId < _maxId));
147   return &_cellIds[_nbDownCells * cellId];
148 }
149
150 /*! get a list of vtk cell types associated to downward entities of a given cell, in the same order
151  * than the downward entities id's list (@see SMDS_Downward::getDownCells).
152  *
153  * @param cellId index of the cell in the downward structure relative to a vtk cell type.
154  * @return table of downward entities types.
155  */
156 const unsigned char* SMDS_Downward::getDownTypes(int cellId)
157 {
158   return &_cellTypes[0];
159 }
160
161 /*! add a downward entity of dimension n-1 (cell or node) to a given cell.
162  * Actual implementation is done in derived methods.
163  * @param cellId index of the parent cell (dimension n) in the downward structure relative to a vtk cell type.
164  * @param lowCellId index of the children cell to add (dimension n-1)
165  * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
166  */
167 void SMDS_Downward::addDownCell(int cellId, int lowCellId, unsigned char aType)
168 {
169   ASSERT(0); // must be re-implemented in derived class
170 }
171
172 /*! add a downward entity of dimension n+1 to a given cell.
173  * Actual implementation is done in derived methods.
174  * @param cellId index of the children cell (dimension n) in the downward structure relative to a vtk cell type.
175  * @param upCellId index of the parent cell to add (dimension n+1)
176  * @param aType vtk cell type of the cell to add (needed to find the SMDS_Downward structure containing the cell to add).
177  */
178 void SMDS_Downward::addUpCell(int cellId, int upCellId, unsigned char aType)
179 {
180   ASSERT(0); // must be re-implemented in derived class
181 }
182
183 int SMDS_Downward::getNodeSet(int cellId, int* nodeSet)
184 {
185   return 0;
186 }
187
188 // ---------------------------------------------------------------------------
189
190 SMDS_Down1D::SMDS_Down1D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
191   SMDS_Downward(grid, nbDownCells)
192 {
193   _upCellIdsVector.clear();
194   _upCellTypesVector.clear();
195   _upCellIds.clear();
196   _upCellTypes.clear();
197   _upCellIndex.clear();
198 }
199
200 SMDS_Down1D::~SMDS_Down1D()
201 {
202 }
203
204 /*! clear vectors used to reference 2D cells containing the edge
205  *
206  * @param cellId
207  */
208 void SMDS_Down1D::initCell(int cellId)
209 {
210   _upCellIdsVector[cellId].clear();
211   _upCellTypesVector[cellId].clear();
212 }
213
214 /*! Resize the downward connectivity storage vector if needed.
215  *
216  * @param nbElems total number of elements of the same type required
217  */
218 void SMDS_Down1D::allocate(int nbElems)
219 {
220   if (nbElems >= _vtkCellIds.size())
221     {
222       _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
223       _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
224       _upCellIdsVector.resize(nbElems + SMDS_Mesh::chunkSize);
225       _upCellTypesVector.resize(nbElems + SMDS_Mesh::chunkSize);
226     }
227 }
228
229 void SMDS_Down1D::compactStorage()
230 {
231   _cellIds.resize(_nbDownCells * _maxId);
232   _vtkCellIds.resize(_maxId);
233
234   int sizeUpCells = 0;
235   for (int i = 0; i < _maxId; i++)
236     sizeUpCells += _upCellIdsVector[i].size();
237   _upCellIds.resize(sizeUpCells, -1);
238   _upCellTypes.resize(sizeUpCells);
239   _upCellIndex.resize(_maxId + 1, -1); // id and types of rank i correspond to [ _upCellIndex[i], _upCellIndex[i+1] [
240   int current = 0;
241   for (int i = 0; i < _maxId; i++)
242     {
243       _upCellIndex[i] = current;
244       for (int j = 0; j < _upCellIdsVector[i].size(); j++)
245         {
246           _upCellIds[current] = _upCellIdsVector[i][j];
247           _upCellTypes[current] = _upCellTypesVector[i][j];
248           current++;
249         }
250     }
251   _upCellIndex[_maxId] = current;
252
253   _upCellIdsVector.clear();
254   _upCellTypesVector.clear();
255 }
256
257 void SMDS_Down1D::addUpCell(int cellId, int upCellId, unsigned char aType)
258 {
259   //ASSERT((cellId >=0) && (cellId < _maxId));
260   int nbFaces = _upCellIdsVector[cellId].size();
261   for (int i = 0; i < nbFaces; i++)
262     {
263       if ((_upCellIdsVector[cellId][i] == upCellId) && (_upCellTypesVector[cellId][i] == aType))
264         {
265           return; // already done
266         }
267     }
268   _upCellIdsVector[cellId].push_back(upCellId);
269   _upCellTypesVector[cellId].push_back(aType);
270 }
271
272 int SMDS_Down1D::getNumberOfUpCells(int cellId)
273 {
274   //ASSERT((cellId >=0) && (cellId < _maxId));
275   return _upCellIndex[cellId + 1] - _upCellIndex[cellId];
276 }
277
278 const int* SMDS_Down1D::getUpCells(int cellId)
279 {
280   //ASSERT((cellId >=0) && (cellId < _maxId));
281   return &_upCellIds[_upCellIndex[cellId]];
282 }
283
284 const unsigned char* SMDS_Down1D::getUpTypes(int cellId)
285 {
286   //ASSERT((cellId >=0) && (cellId < _maxId));
287   return &_upCellTypes[_upCellIndex[cellId]];
288 }
289
290 void SMDS_Down1D::getNodeIds(int cellId, std::set<int>& nodeSet)
291 {
292   for (int i = 0; i < _nbDownCells; i++)
293     nodeSet.insert(_cellIds[_nbDownCells * cellId + i]);
294 }
295
296 int SMDS_Down1D::getNodeSet(int cellId, int* nodeSet)
297 {
298   for (int i = 0; i < _nbDownCells; i++)
299     nodeSet[i] = _cellIds[_nbDownCells * cellId + i];
300   return _nbDownCells;
301 }
302
303 void SMDS_Down1D::setNodes(int cellId, int vtkId)
304 {
305   vtkIdType npts = 0;
306   vtkIdType *pts; // will refer to the point id's of the face
307   _grid->GetCellPoints(vtkId, npts, pts);
308   // MESSAGE(vtkId << " " << npts << "  " << _nbDownCells);
309   //ASSERT(npts == _nbDownCells);
310   for (int i = 0; i < npts; i++)
311     {
312       _cellIds[_nbDownCells * cellId + i] = pts[i];
313     }
314 }
315
316 void SMDS_Down1D::setNodes(int cellId, const int* nodeIds)
317 {
318   //ASSERT(nodeIds.size() == _nbDownCells);
319   for (int i = 0; i < _nbDownCells; i++)
320     {
321       _cellIds[_nbDownCells * cellId + i] = nodeIds[i];
322     }
323 }
324
325 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
326  * We keep in the list the cells that contains all the nodes, we keep only volumes and faces.
327  * @param cellId id of the edge in the downward structure
328  * @param vtkIds vector of vtk id's
329  * @return number of vtk cells (size of vector)
330  */
331 int SMDS_Down1D::computeVtkCells(int cellId, std::vector<int>& vtkIds)
332 {
333   vtkIds.clear();
334
335   // --- find all the cells the points belong to, and how many of the points belong to a given cell
336
337   int *pts = &_cellIds[_nbDownCells * cellId];
338   int ncells = this->computeVtkCells(pts, vtkIds);
339   return ncells;
340 }
341
342 /*! Build the list of vtkUnstructuredGrid cells containing the edge.
343  *
344  * @param pts list of points id's defining an edge
345  * @param vtkIds vector of vtk id's
346  * @return number of vtk cells (size of vector)
347  */
348 int SMDS_Down1D::computeVtkCells(int *pts, std::vector<int>& vtkIds)
349 {
350
351   // --- find all the cells the points belong to, and how many of the points belong to a given cell
352
353   int cellIds[1000];
354   int cellCnt[1000];
355   int cnt = 0;
356   for (int i = 0; i < _nbDownCells; i++)
357     {
358       vtkIdType point = pts[i];
359       int numCells = _grid->GetLinks()->GetNcells(point);
360       vtkIdType *cells = _grid->GetLinks()->GetCells(point);
361       for (int j = 0; j < numCells; j++)
362         {
363           int vtkCellId = cells[j];
364           bool found = false;
365           for (int k = 0; k < cnt; k++)
366             {
367               if (cellIds[k] == vtkCellId)
368                 {
369                   cellCnt[k] += 1;
370                   found = true;
371                   break;
372                 }
373             }
374           if (!found)
375             {
376               cellIds[cnt] = vtkCellId;
377               cellCnt[cnt] = 1;
378               // TODO ASSERT(cnt<1000);
379               cnt++;
380             }
381         }
382     }
383
384   // --- find the face and volume cells: they contains all the points and are of type volume or face
385
386   int ncells = 0;
387   for (int i = 0; i < cnt; i++)
388     {
389       if (cellCnt[i] == _nbDownCells)
390         {
391           int vtkElemId = cellIds[i];
392           int vtkType = _grid->GetCellType(vtkElemId);
393           if (SMDS_Downward::getCellDimension(vtkType) > 1)
394             {
395               vtkIds.push_back(vtkElemId);
396               ncells++;
397             }
398         }
399     }
400
401   return ncells;
402 }
403
404 /*! Build the list of downward faces from a list of vtk cells.
405  *
406  * @param cellId id of the edge in the downward structure
407  * @param vtkIds vector of vtk id's
408  * @param downFaces vector of face id's in downward structures
409  * @param downTypes vector of face types
410  * @return number of downward faces
411  */
412 int SMDS_Down1D::computeFaces(int cellId, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
413 {
414   int *pts = &_cellIds[_nbDownCells * cellId];
415   int nbFaces = this->computeFaces(pts, vtkIds, nbcells, downFaces, downTypes);
416   return nbFaces;
417 }
418
419 /*! Build the list of downward faces from a list of vtk cells.
420  *
421  * @param pts list of points id's defining an edge
422  * @param vtkIds vector of vtk id's
423  * @param downFaces vector of face id's in downward structures
424  * @param downTypes vector of face types
425  * @return number of downward faces
426  */
427 int SMDS_Down1D::computeFaces(int* pts, int* vtkIds, int nbcells, int* downFaces, unsigned char* downTypes)
428 {
429   int cnt = 0;
430   for (int i = 0; i < nbcells; i++)
431     {
432       int vtkId = vtkIds[i];
433       int vtkType = _grid->GetCellType(vtkId);
434       if (SMDS_Downward::getCellDimension(vtkType) == 2)
435         {
436           int faceId = _grid->CellIdToDownId(vtkId);
437           downFaces[cnt] = faceId;
438           downTypes[cnt] = vtkType;
439           cnt++;
440         }
441       else if (SMDS_Downward::getCellDimension(vtkType) == 3)
442         {
443           int volId = _grid->CellIdToDownId(vtkId);
444           SMDS_Downward * downvol = _grid->getDownArray(vtkType);
445           //const int *downIds = downvol->getDownCells(volId);
446           const unsigned char* downTypesVol = downvol->getDownTypes(volId);
447           int nbFaces = downvol->getNumberOfDownCells(volId);
448           const int* faceIds = downvol->getDownCells(volId);
449           for (int n = 0; n < nbFaces; n++)
450             {
451               SMDS_Down2D *downFace = static_cast<SMDS_Down2D*> (_grid->getDownArray(downTypesVol[n]));
452               bool isInFace = downFace->isInFace(faceIds[n], pts, _nbDownCells);
453               if (isInFace)
454                 {
455                   bool alreadySet = false;
456                   for (int k = 0; k < cnt; k++)
457                     if (faceIds[n] == downFaces[k])
458                       {
459                         alreadySet = true;
460                         break;
461                       }
462                   if (!alreadySet)
463                     {
464                       downFaces[cnt] = faceIds[n];
465                       downTypes[cnt] = downTypesVol[n];
466                       cnt++;
467                     }
468                 }
469             }
470         }
471     }
472   return cnt;
473 }
474
475 // ---------------------------------------------------------------------------
476
477 SMDS_Down2D::SMDS_Down2D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
478   SMDS_Downward(grid, nbDownCells)
479 {
480   _upCellIds.clear();
481   _upCellTypes.clear();
482   _tempNodes.clear();
483   _nbNodes = 0;
484 }
485
486 SMDS_Down2D::~SMDS_Down2D()
487 {
488 }
489
490 int SMDS_Down2D::getNumberOfUpCells(int cellId)
491 {
492   int nbup = 0;
493   if (_upCellIds[2 * cellId] >= 0)
494     nbup++;
495   if (_upCellIds[2 * cellId + 1] >= 0)
496     nbup++;
497   return nbup;
498 }
499
500 const int* SMDS_Down2D::getUpCells(int cellId)
501 {
502   //ASSERT((cellId >=0) && (cellId < _maxId));
503   return &_upCellIds[2 * cellId];
504 }
505
506 const unsigned char* SMDS_Down2D::getUpTypes(int cellId)
507 {
508   //ASSERT((cellId >=0) && (cellId < _maxId));
509   return &_upCellTypes[2 * cellId];
510 }
511
512 void SMDS_Down2D::getNodeIds(int cellId, std::set<int>& nodeSet)
513 {
514   for (int i = 0; i < _nbDownCells; i++)
515     {
516       int downCellId = _cellIds[_nbDownCells * cellId + i];
517       unsigned char cellType = _cellTypes[i];
518       this->_grid->getDownArray(cellType)->getNodeIds(downCellId, nodeSet);
519     }
520 }
521
522 /*! Find in vtkUnstructuredGrid the volumes containing a face already stored in vtkUnstructuredGrid.
523  * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
524  * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
525  * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
526  * @param cellId the face cell id in vkUnstructuredGrid
527  * @param ids a couple of vtkId, initialized at -1 (no parent volume)
528  * @return number of volumes (0, 1 or 2)
529  */
530 int SMDS_Down2D::computeVolumeIds(int cellId, int* ids)
531 {
532   // --- find point id's of the face
533
534   vtkIdType npts = 0;
535   vtkIdType *pts; // will refer to the point id's of the face
536   _grid->GetCellPoints(cellId, npts, pts);
537   vector<int> nodes;
538   for (int i = 0; i < npts; i++)
539     nodes.push_back(pts[i]);
540   int nvol = this->computeVolumeIdsFromNodesFace(&nodes[0], npts, ids);
541   return nvol;
542 }
543
544 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
545  * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
546  * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
547  * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
548  * @param faceByNodes
549  * @param ids a couple of vtkId, initialized at -1 (no parent volume)
550  * @return number of volumes (0, 1 or 2)
551  */
552 int SMDS_Down2D::computeVolumeIds(ElemByNodesType& faceByNodes, int* ids)
553 {
554   int nvol = this->computeVolumeIdsFromNodesFace(&faceByNodes.nodeIds[0], faceByNodes.nbNodes, ids);
555   return nvol;
556 }
557
558 /*! Find in vtkUnstructuredGrid the volumes containing a face described by it's nodes
559  * Search the volumes containing a face, to store the info in SMDS_Down2D for later uses
560  * with SMDS_Down2D::getUpCells and SMDS_Down2D::getUpTypes.
561  * A face belongs to 0, 1 or 2 volumes, identified by their id in vtkUnstructuredGrid.
562  * @param pts array of vtk node id's
563  * @param npts number of nodes
564  * @param ids
565  * @return number of volumes (0, 1 or 2)
566  */
567 int SMDS_Down2D::computeVolumeIdsFromNodesFace(int* pts, int npts, int* ids)
568 {
569
570   // --- find all the cells the points belong to, and how many of the points belong to a given cell
571
572   int cellIds[1000];
573   int cellCnt[1000];
574   int cnt = 0;
575   for (int i = 0; i < npts; i++)
576     {
577       vtkIdType point = pts[i];
578       int numCells = _grid->GetLinks()->GetNcells(point);
579       //MESSAGE("cells pour " << i << " " << numCells);
580       vtkIdType *cells = _grid->GetLinks()->GetCells(point);
581       for (int j = 0; j < numCells; j++)
582         {
583           int vtkCellId = cells[j];
584           bool found = false;
585           for (int k = 0; k < cnt; k++)
586             {
587               if (cellIds[k] == vtkCellId)
588                 {
589                   cellCnt[k] += 1;
590                   found = true;
591                   break;
592                 }
593             }
594           if (!found)
595             {
596               cellIds[cnt] = vtkCellId;
597               cellCnt[cnt] = 1;
598               // TODO ASSERT(cnt<1000);
599               cnt++;
600             }
601         }
602     }
603
604   // --- find the volume cells: they contains all the points and are of type volume
605
606   int nvol = 0;
607   for (int i = 0; i < cnt; i++)
608     {
609       //MESSAGE("cell " << cellIds[i] << " points " << cellCnt[i]);
610       if (cellCnt[i] == npts)
611         {
612           int vtkElemId = cellIds[i];
613           int vtkType = _grid->GetCellType(vtkElemId);
614           if (SMDS_Downward::getCellDimension(vtkType) == 3)
615             {
616               ids[nvol] = vtkElemId; // store the volume id in given vector
617               nvol++;
618             }
619         }
620       if (nvol == 2)
621         break;
622     }
623
624   return nvol;
625 }
626
627 void SMDS_Down2D::setTempNodes(int cellId, int vtkId)
628 {
629   vtkIdType npts = 0;
630   vtkIdType *pts; // will refer to the point id's of the face
631   _grid->GetCellPoints(vtkId, npts, pts);
632   // MESSAGE(vtkId << " " << npts << "  " << _nbNodes);
633   //ASSERT(npts == _nbNodes);
634   for (int i = 0; i < npts; i++)
635     {
636       _tempNodes[_nbNodes * cellId + i] = pts[i];
637     }
638 }
639
640 void SMDS_Down2D::setTempNodes(int cellId, ElemByNodesType& faceByNodes)
641 {
642   for (int i = 0; i < faceByNodes.nbNodes; i++)
643     _tempNodes[_nbNodes * cellId + i] = faceByNodes.nodeIds[i];
644 }
645
646 /*! Find if all the nodes belongs to the face.
647  *
648  * @param cellId the face cell Id
649  * @param nodeSet set of node id's to be found in the face list of nodes
650  * @return
651  */
652 bool SMDS_Down2D::isInFace(int cellId, int *pts, int npts)
653 {
654   int nbFound = 0;
655   int *nodes = &_tempNodes[_nbNodes * cellId];
656   for (int j = 0; j < npts; j++)
657     {
658       int point = pts[j];
659       for (int i = 0; i < _nbNodes; i++)
660         {
661           if (nodes[i] == point)
662             {
663               nbFound++;
664               break;
665             }
666         }
667     }
668   return (nbFound == npts);
669 }
670
671 /*! Resize the downward connectivity storage vector if needed.
672  *
673  * @param nbElems total number of elements of the same type required
674  */
675 void SMDS_Down2D::allocate(int nbElems)
676 {
677   if (nbElems >= _vtkCellIds.size())
678     {
679       _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
680       _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
681       _upCellIds.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
682       _upCellTypes.resize(2 * (nbElems + SMDS_Mesh::chunkSize), -1);
683       _tempNodes.resize(_nbNodes * (nbElems + SMDS_Mesh::chunkSize), -1);
684     }
685 }
686
687 void SMDS_Down2D::compactStorage()
688 {
689   _cellIds.resize(_nbDownCells * _maxId);
690   _upCellIds.resize(2 * _maxId);
691   _upCellTypes.resize(2 * _maxId);
692   _vtkCellIds.resize(_maxId);
693   _tempNodes.clear();
694 }
695
696 void SMDS_Down2D::addUpCell(int cellId, int upCellId, unsigned char aType)
697 {
698   //ASSERT((cellId >=0)&& (cellId < _maxId));
699   int *vols = &_upCellIds[2 * cellId];
700   unsigned char *types = &_upCellTypes[2 * cellId];
701   for (int i = 0; i < 2; i++)
702     {
703       if (vols[i] < 0)
704         {
705           vols[i] = upCellId; // use non affected volume
706           types[i] = aType;
707           return;
708         }
709       if ((vols[i] == upCellId) && (types[i] == aType)) // already done
710         return;
711     }
712   ASSERT(0);
713 }
714
715 int SMDS_Down2D::getNodeSet(int cellId, int* nodeSet)
716 {
717   for (int i = 0; i < _nbNodes; i++)
718     nodeSet[i] = _tempNodes[_nbNodes * cellId + i];
719   return _nbNodes;
720 }
721
722 int SMDS_Down2D::FindEdgeByNodes(int cellId, ElemByNodesType& edgeByNodes)
723 {
724   int *edges = &_cellIds[_nbDownCells * cellId];
725   for (int i = 0; i < _nbDownCells; i++)
726     {
727       if ((edges[i] >= 0) && (edgeByNodes.vtkType == _cellTypes[i]))
728         {
729           int nodeSet[3];
730           int npts = this->_grid->getDownArray(edgeByNodes.vtkType)->getNodeSet(edges[i], nodeSet);
731           bool found = false;
732           for (int j = 0; j < npts; j++)
733             {
734               int point = edgeByNodes.nodeIds[j];
735               found = false;
736               for (int k = 0; k < npts; k++)
737                 {
738                   if (nodeSet[k] == point)
739                     {
740                       found = true;
741                       break;
742                     }
743                 }
744               if (!found)
745                 break;
746             }
747           if (found)
748             return edges[i];
749         }
750     }
751   return -1;
752 }
753
754 // ---------------------------------------------------------------------------
755
756 SMDS_Down3D::SMDS_Down3D(SMDS_UnstructuredGrid *grid, int nbDownCells) :
757   SMDS_Downward(grid, nbDownCells)
758 {
759 }
760
761 SMDS_Down3D::~SMDS_Down3D()
762 {
763 }
764
765 void SMDS_Down3D::allocate(int nbElems)
766 {
767   if (nbElems >= _vtkCellIds.size())
768     {
769       _cellIds.resize(_nbDownCells * (nbElems + SMDS_Mesh::chunkSize), -1);
770       _vtkCellIds.resize(nbElems + SMDS_Mesh::chunkSize, -1);
771     }
772 }
773
774 void SMDS_Down3D::compactStorage()
775 {
776   // nothing to do, size was known before
777 }
778
779 int SMDS_Down3D::getNumberOfUpCells(int cellId)
780 {
781   return 0;
782 }
783
784 const int* SMDS_Down3D::getUpCells(int cellId)
785 {
786   return 0;
787 }
788
789 const unsigned char* SMDS_Down3D::getUpTypes(int cellId)
790 {
791   return 0;
792 }
793
794 void SMDS_Down3D::getNodeIds(int cellId, std::set<int>& nodeSet)
795 {
796   int vtkId = this->_vtkCellIds[cellId];
797   vtkIdType npts = 0;
798   vtkIdType *nodes; // will refer to the point id's of the volume
799   _grid->GetCellPoints(vtkId, npts, nodes);
800   for (int i = 0; i < npts; i++)
801     nodeSet.insert(nodes[i]);
802 }
803
804 int SMDS_Down3D::FindFaceByNodes(int cellId, ElemByNodesType& faceByNodes)
805 {
806   int *faces = &_cellIds[_nbDownCells * cellId];
807   int npoints = 0;
808
809   for (int i = 0; i < _nbDownCells; i++)
810     {
811       if ((faces[i] >= 0) && (faceByNodes.vtkType == _cellTypes[i]))
812         {
813           if (npoints == 0)
814             npoints = faceByNodes.nbNodes;
815
816           int nodeSet[10];
817           int npts = this->_grid->getDownArray(faceByNodes.vtkType)->getNodeSet(faces[i], nodeSet);
818           if (npts != npoints)
819             continue; // skip this face
820           bool found = false;
821           for (int j = 0; j < npts; j++)
822             {
823               int point = faceByNodes.nodeIds[j];
824               found = false;
825               for (int k = 0; k < npts; k++)
826                 {
827                   if (nodeSet[k] == point)
828                     {
829                       found = true;
830                       break; // point j is in the 2 faces, skip remaining k values
831                     }
832                 }
833               if (!found)
834                 break; // point j is not in the 2 faces, skip the remaining tests
835             }
836           if (found)
837             return faces[i];
838         }
839     }
840   return -1;
841 }
842
843 // ---------------------------------------------------------------------------
844
845 SMDS_DownEdge::SMDS_DownEdge(SMDS_UnstructuredGrid *grid) :
846   SMDS_Down1D(grid, 2)
847 {
848   _cellTypes.push_back(VTK_VERTEX);
849   _cellTypes.push_back(VTK_VERTEX);
850 }
851
852 SMDS_DownEdge::~SMDS_DownEdge()
853 {
854 }
855
856 // ---------------------------------------------------------------------------
857
858 SMDS_DownQuadEdge::SMDS_DownQuadEdge(SMDS_UnstructuredGrid *grid) :
859   SMDS_Down1D(grid, 3)
860 {
861   _cellTypes.push_back(VTK_VERTEX);
862   _cellTypes.push_back(VTK_VERTEX);
863   _cellTypes.push_back(VTK_VERTEX);
864 }
865
866 SMDS_DownQuadEdge::~SMDS_DownQuadEdge()
867 {
868 }
869
870 // ---------------------------------------------------------------------------
871
872 SMDS_DownTriangle::SMDS_DownTriangle(SMDS_UnstructuredGrid *grid) :
873   SMDS_Down2D(grid, 3)
874 {
875   _cellTypes.push_back(VTK_LINE);
876   _cellTypes.push_back(VTK_LINE);
877   _cellTypes.push_back(VTK_LINE);
878   _nbNodes = 3;
879 }
880
881 SMDS_DownTriangle::~SMDS_DownTriangle()
882 {
883 }
884
885 void SMDS_DownTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
886 {
887   int *nodes = &_tempNodes[_nbNodes * cellId];
888   edgesWithNodes.nbElems = 3;
889
890   edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
891   edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
892   edgesWithNodes.elems[0].nbNodes = 2;
893   edgesWithNodes.elems[0].vtkType = VTK_LINE;
894
895   edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
896   edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
897   edgesWithNodes.elems[1].nbNodes = 2;
898   edgesWithNodes.elems[1].vtkType = VTK_LINE;
899
900   edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
901   edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
902   edgesWithNodes.elems[2].nbNodes = 2;
903   edgesWithNodes.elems[2].vtkType = VTK_LINE;
904 }
905
906 void SMDS_DownTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
907 {
908   //ASSERT((cellId >=0)&& (cellId < _maxId));
909   //ASSERT(aType == VTK_LINE);
910   int *faces = &_cellIds[_nbDownCells * cellId];
911   for (int i = 0; i < _nbDownCells; i++)
912     {
913       if (faces[i] < 0)
914         {
915           faces[i] = lowCellId;
916           return;
917         }
918       if (faces[i] == lowCellId)
919         return;
920     }
921   ASSERT(0);
922 }
923
924 // ---------------------------------------------------------------------------
925
926 SMDS_DownQuadTriangle::SMDS_DownQuadTriangle(SMDS_UnstructuredGrid *grid) :
927   SMDS_Down2D(grid, 3)
928 {
929   _cellTypes.push_back(VTK_QUADRATIC_EDGE);
930   _cellTypes.push_back(VTK_QUADRATIC_EDGE);
931   _cellTypes.push_back(VTK_QUADRATIC_EDGE);
932   _nbNodes = 6;
933 }
934
935 SMDS_DownQuadTriangle::~SMDS_DownQuadTriangle()
936 {
937 }
938
939 void SMDS_DownQuadTriangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
940 {
941   int *nodes = &_tempNodes[_nbNodes * cellId];
942   edgesWithNodes.nbElems = 3;
943
944   edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
945   edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
946   edgesWithNodes.elems[0].nodeIds[2] = nodes[3];
947   edgesWithNodes.elems[0].nbNodes = 3;
948   edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
949
950   edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
951   edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
952   edgesWithNodes.elems[1].nodeIds[2] = nodes[4];
953   edgesWithNodes.elems[1].nbNodes = 3;
954   edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
955
956   edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
957   edgesWithNodes.elems[2].nodeIds[1] = nodes[0];
958   edgesWithNodes.elems[2].nodeIds[2] = nodes[5];
959   edgesWithNodes.elems[2].nbNodes = 3;
960   edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
961 }
962
963 void SMDS_DownQuadTriangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
964 {
965   //ASSERT((cellId >=0)&& (cellId < _maxId));
966   //ASSERT(aType == VTK_QUADRATIC_EDGE);
967   int *edges = &_cellIds[_nbDownCells * cellId];
968   for (int i = 0; i < _nbDownCells; i++)
969     {
970       if (edges[i] < 0)
971         {
972           edges[i] = lowCellId;
973           return;
974         }
975       if (edges[i] == lowCellId)
976         return;
977     }
978   ASSERT(0);
979 }
980
981 // ---------------------------------------------------------------------------
982
983 SMDS_DownQuadrangle::SMDS_DownQuadrangle(SMDS_UnstructuredGrid *grid) :
984   SMDS_Down2D(grid, 4)
985 {
986   _cellTypes.push_back(VTK_LINE);
987   _cellTypes.push_back(VTK_LINE);
988   _cellTypes.push_back(VTK_LINE);
989   _cellTypes.push_back(VTK_LINE);
990   _nbNodes = 4;
991 }
992
993 SMDS_DownQuadrangle::~SMDS_DownQuadrangle()
994 {
995 }
996
997 void SMDS_DownQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
998 {
999   int *nodes = &_tempNodes[_nbNodes * cellId];
1000   edgesWithNodes.nbElems = 4;
1001
1002   edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1003   edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1004   edgesWithNodes.elems[0].nbNodes = 2;
1005   edgesWithNodes.elems[0].vtkType = VTK_LINE;
1006
1007   edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1008   edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1009   edgesWithNodes.elems[1].nbNodes = 2;
1010   edgesWithNodes.elems[1].vtkType = VTK_LINE;
1011
1012   edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1013   edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1014   edgesWithNodes.elems[2].nbNodes = 2;
1015   edgesWithNodes.elems[2].vtkType = VTK_LINE;
1016
1017   edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1018   edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1019   edgesWithNodes.elems[3].nbNodes = 2;
1020   edgesWithNodes.elems[3].vtkType = VTK_LINE;
1021 }
1022
1023 void SMDS_DownQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1024 {
1025   //ASSERT((cellId >=0)&& (cellId < _maxId));
1026   //ASSERT(aType == VTK_LINE);
1027   int *faces = &_cellIds[_nbDownCells * cellId];
1028   for (int i = 0; i < _nbDownCells; i++)
1029     {
1030       if (faces[i] < 0)
1031         {
1032           faces[i] = lowCellId;
1033           return;
1034         }
1035       if (faces[i] == lowCellId)
1036         return;
1037     }
1038   ASSERT(0);
1039 }
1040
1041 // ---------------------------------------------------------------------------
1042
1043 SMDS_DownQuadQuadrangle::SMDS_DownQuadQuadrangle(SMDS_UnstructuredGrid *grid) :
1044   SMDS_Down2D(grid, 4)
1045 {
1046   _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1047   _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1048   _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1049   _cellTypes.push_back(VTK_QUADRATIC_EDGE);
1050   _nbNodes = 8;
1051 }
1052
1053 SMDS_DownQuadQuadrangle::~SMDS_DownQuadQuadrangle()
1054 {
1055 }
1056
1057 void SMDS_DownQuadQuadrangle::computeEdgesWithNodes(int cellId, ListElemByNodesType& edgesWithNodes)
1058 {
1059   int *nodes = &_tempNodes[_nbNodes * cellId];
1060   edgesWithNodes.nbElems = 4;
1061
1062   edgesWithNodes.elems[0].nodeIds[0] = nodes[0];
1063   edgesWithNodes.elems[0].nodeIds[1] = nodes[1];
1064   edgesWithNodes.elems[0].nodeIds[2] = nodes[4];
1065   edgesWithNodes.elems[0].nbNodes = 3;
1066   edgesWithNodes.elems[0].vtkType = VTK_QUADRATIC_EDGE;
1067
1068   edgesWithNodes.elems[1].nodeIds[0] = nodes[1];
1069   edgesWithNodes.elems[1].nodeIds[1] = nodes[2];
1070   edgesWithNodes.elems[1].nodeIds[2] = nodes[5];
1071   edgesWithNodes.elems[1].nbNodes = 3;
1072   edgesWithNodes.elems[1].vtkType = VTK_QUADRATIC_EDGE;
1073
1074   edgesWithNodes.elems[2].nodeIds[0] = nodes[2];
1075   edgesWithNodes.elems[2].nodeIds[1] = nodes[3];
1076   edgesWithNodes.elems[2].nodeIds[2] = nodes[6];
1077   edgesWithNodes.elems[2].nbNodes = 3;
1078   edgesWithNodes.elems[2].vtkType = VTK_QUADRATIC_EDGE;
1079
1080   edgesWithNodes.elems[3].nodeIds[0] = nodes[3];
1081   edgesWithNodes.elems[3].nodeIds[1] = nodes[0];
1082   edgesWithNodes.elems[3].nodeIds[2] = nodes[7];
1083   edgesWithNodes.elems[3].nbNodes = 3;
1084   edgesWithNodes.elems[3].vtkType = VTK_QUADRATIC_EDGE;
1085 }
1086
1087 void SMDS_DownQuadQuadrangle::addDownCell(int cellId, int lowCellId, unsigned char aType)
1088 {
1089   //ASSERT((cellId >=0)&& (cellId < _maxId));
1090   //ASSERT(aType == VTK_QUADRATIC_EDGE);
1091   int *faces = &_cellIds[_nbDownCells * cellId];
1092   for (int i = 0; i < _nbDownCells; i++)
1093     {
1094       if (faces[i] < 0)
1095         {
1096           faces[i] = lowCellId;
1097           return;
1098         }
1099       if (faces[i] == lowCellId)
1100         return;
1101     }
1102   ASSERT(0);
1103 }
1104
1105 // ---------------------------------------------------------------------------
1106
1107 SMDS_DownTetra::SMDS_DownTetra(SMDS_UnstructuredGrid *grid) :
1108   SMDS_Down3D(grid, 4)
1109 {
1110   _cellTypes.push_back(VTK_TRIANGLE);
1111   _cellTypes.push_back(VTK_TRIANGLE);
1112   _cellTypes.push_back(VTK_TRIANGLE);
1113   _cellTypes.push_back(VTK_TRIANGLE);
1114 }
1115
1116 SMDS_DownTetra::~SMDS_DownTetra()
1117 {
1118 }
1119
1120 void SMDS_DownTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1121 {
1122   set<int> setNodes;
1123   setNodes.clear();
1124   for (int i = 0; i < orderedNodes.size(); i++)
1125     setNodes.insert(orderedNodes[i]);
1126   //MESSAGE("cellId = " << cellId);
1127
1128   vtkIdType npts = 0;
1129   vtkIdType *nodes; // will refer to the point id's of the volume
1130   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1131
1132   set<int> tofind;
1133   int ids[12] = { 0, 1, 2,  0, 3, 1,  2, 3, 0,   1, 3, 2 };
1134 //int ids[12] = { 2, 1, 0,  1, 3, 0,  0, 3, 2,   2, 3, 1 };
1135   for (int k = 0; k < 4; k++)
1136     {
1137       tofind.clear();
1138       for (int i = 0; i < 3; i++)
1139         tofind.insert(nodes[ids[3 * k + i]]);
1140       if (setNodes == tofind)
1141         {
1142           for (int i = 0; i < 3; i++)
1143             orderedNodes[i] = nodes[ids[3 * k + i]];
1144           return;
1145         }
1146     }
1147   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1148   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1149   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1150 }
1151
1152 void SMDS_DownTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1153 {
1154   //ASSERT((cellId >=0)&& (cellId < _maxId));
1155   //ASSERT(aType == VTK_TRIANGLE);
1156   int *faces = &_cellIds[_nbDownCells * cellId];
1157   for (int i = 0; i < _nbDownCells; i++)
1158     {
1159       if (faces[i] < 0)
1160         {
1161           faces[i] = lowCellId;
1162           return;
1163         }
1164       if (faces[i] == lowCellId)
1165         return;
1166     }
1167   ASSERT(0);
1168 }
1169
1170 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's
1171  * The linear tetrahedron is defined by four points.
1172  * @see vtkTetra.h in Filtering.
1173  * @param cellId volumeId in vtkUnstructuredGrid
1174  * @param facesWithNodes vector of face descriptors to be filled
1175  */
1176 void SMDS_DownTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1177 {
1178   // --- find point id's of the volume
1179
1180   vtkIdType npts = 0;
1181   vtkIdType *nodes; // will refer to the point id's of the volume
1182   _grid->GetCellPoints(cellId, npts, nodes);
1183
1184   // --- create all the ordered list of node id's for each face
1185
1186   facesWithNodes.nbElems = 4;
1187
1188   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1189   facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1190   facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1191   facesWithNodes.elems[0].nbNodes = 3;
1192   facesWithNodes.elems[0].vtkType = VTK_TRIANGLE;
1193
1194   facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1195   facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1196   facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1197   facesWithNodes.elems[1].nbNodes = 3;
1198   facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1199
1200   facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1201   facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1202   facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1203   facesWithNodes.elems[2].nbNodes = 3;
1204   facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1205
1206   facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1207   facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1208   facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1209   facesWithNodes.elems[3].nbNodes = 3;
1210   facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1211 }
1212
1213 // ---------------------------------------------------------------------------
1214
1215 SMDS_DownQuadTetra::SMDS_DownQuadTetra(SMDS_UnstructuredGrid *grid) :
1216   SMDS_Down3D(grid, 4)
1217 {
1218   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1219   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1220   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1221   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1222 }
1223
1224 SMDS_DownQuadTetra::~SMDS_DownQuadTetra()
1225 {
1226 }
1227
1228 void SMDS_DownQuadTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1229 {
1230   set<int> setNodes;
1231   setNodes.clear();
1232   for (int i = 0; i < orderedNodes.size(); i++)
1233     setNodes.insert(orderedNodes[i]);
1234   //MESSAGE("cellId = " << cellId);
1235
1236   vtkIdType npts = 0;
1237   vtkIdType *nodes; // will refer to the point id's of the volume
1238   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1239
1240   set<int> tofind;
1241   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 };
1242 //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 };
1243   for (int k = 0; k < 4; k++)
1244     {
1245       tofind.clear();
1246       for (int i = 0; i < 6; i++)
1247         tofind.insert(nodes[ids[6 * k + i]]);
1248       if (setNodes == tofind)
1249         {
1250           for (int i = 0; i < 6; i++)
1251             orderedNodes[i] = nodes[ids[6 * k + i]];
1252           return;
1253         }
1254     }
1255   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1256   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1257   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1258 }
1259
1260 void SMDS_DownQuadTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
1261 {
1262   //ASSERT((cellId >=0)&& (cellId < _maxId));
1263   //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1264   int *faces = &_cellIds[_nbDownCells * cellId];
1265   for (int i = 0; i < _nbDownCells; i++)
1266     {
1267       if (faces[i] < 0)
1268         {
1269           faces[i] = lowCellId;
1270           return;
1271         }
1272       if (faces[i] == lowCellId)
1273         return;
1274     }
1275   ASSERT(0);
1276 }
1277
1278 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's
1279  * The ordering of the ten points defining the quadratic tetrahedron cell is point id's (0-3,4-9)
1280  * where id's 0-3 are the four tetrahedron vertices;
1281  * 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).
1282  * @see vtkQuadraticTetra.h in Filtering.
1283  * @param cellId volumeId in vtkUnstructuredGrid
1284  * @param facesWithNodes vector of face descriptors to be filled
1285  */
1286 void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1287 {
1288   // --- find point id's of the volume
1289
1290   vtkIdType npts = 0;
1291   vtkIdType *nodes; // will refer to the point id's of the volume
1292   _grid->GetCellPoints(cellId, npts, nodes);
1293
1294   // --- create all the ordered list of node id's for each face
1295
1296   facesWithNodes.nbElems = 4;
1297
1298   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1299   facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1300   facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1301   facesWithNodes.elems[0].nodeIds[3] = nodes[4];
1302   facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1303   facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1304   facesWithNodes.elems[0].nbNodes = 6;
1305   facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_TRIANGLE;
1306
1307   facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1308   facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1309   facesWithNodes.elems[1].nodeIds[2] = nodes[3];
1310   facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1311   facesWithNodes.elems[1].nodeIds[4] = nodes[8];
1312   facesWithNodes.elems[1].nodeIds[5] = nodes[7];
1313   facesWithNodes.elems[1].nbNodes = 6;
1314   facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1315
1316   facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1317   facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1318   facesWithNodes.elems[2].nodeIds[2] = nodes[3];
1319   facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1320   facesWithNodes.elems[2].nodeIds[4] = nodes[9];
1321   facesWithNodes.elems[2].nodeIds[5] = nodes[7];
1322   facesWithNodes.elems[2].nbNodes = 6;
1323   facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1324
1325   facesWithNodes.elems[3].nodeIds[0] = nodes[1];
1326   facesWithNodes.elems[3].nodeIds[1] = nodes[2];
1327   facesWithNodes.elems[3].nodeIds[2] = nodes[3];
1328   facesWithNodes.elems[3].nodeIds[3] = nodes[5];
1329   facesWithNodes.elems[3].nodeIds[4] = nodes[9];
1330   facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1331   facesWithNodes.elems[3].nbNodes = 6;
1332   facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1333 }
1334
1335 // ---------------------------------------------------------------------------
1336
1337 SMDS_DownPyramid::SMDS_DownPyramid(SMDS_UnstructuredGrid *grid) :
1338   SMDS_Down3D(grid, 5)
1339 {
1340   _cellTypes.push_back(VTK_QUAD);
1341   _cellTypes.push_back(VTK_TRIANGLE);
1342   _cellTypes.push_back(VTK_TRIANGLE);
1343   _cellTypes.push_back(VTK_TRIANGLE);
1344   _cellTypes.push_back(VTK_TRIANGLE);
1345 }
1346
1347 SMDS_DownPyramid::~SMDS_DownPyramid()
1348 {
1349 }
1350
1351 void SMDS_DownPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1352 {
1353   set<int> setNodes;
1354   setNodes.clear();
1355   for (int i = 0; i < orderedNodes.size(); i++)
1356     setNodes.insert(orderedNodes[i]);
1357   //MESSAGE("cellId = " << cellId);
1358
1359   vtkIdType npts = 0;
1360   vtkIdType *nodes; // will refer to the point id's of the volume
1361   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1362
1363   set<int> tofind;
1364   int ids[16] = { 0, 1, 2, 3,   0, 3, 4,   3, 2, 4,   2, 1, 4,   1, 0, 4 };
1365
1366   // Quadrangular face
1367   tofind.clear();
1368   for (int i = 0; i < 4; i++)
1369     tofind.insert(nodes[ids[i]]);
1370   if (setNodes == tofind)
1371     {
1372       for (int i = 0; i < 4; i++)
1373         orderedNodes[i] = nodes[ids[i]];
1374       return;
1375     }
1376   // Triangular faces
1377   for (int k = 0; k < 4; k++)
1378     {
1379       tofind.clear();
1380       for (int i = 0; i < 3; i++)
1381         tofind.insert(nodes[ids[4 + 3 * k + i]]);
1382       if (setNodes == tofind)
1383         {
1384           for (int i = 0; i < 3; i++)
1385             orderedNodes[i] = nodes[ids[4 + 3 * k + i]];
1386           return;
1387         }
1388     }
1389   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1390   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1391   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1392 }
1393
1394 void SMDS_DownPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1395 {
1396   //ASSERT((cellId >=0) && (cellId < _maxId));
1397   int *faces = &_cellIds[_nbDownCells * cellId];
1398   if (aType == VTK_QUAD)
1399     {
1400       if (faces[0] < 0)
1401         {
1402           faces[0] = lowCellId;
1403           return;
1404         }
1405       if (faces[0] == lowCellId)
1406         return;
1407     }
1408   else
1409     {
1410       //ASSERT(aType == VTK_TRIANGLE);
1411       for (int i = 1; i < _nbDownCells; i++)
1412         {
1413           if (faces[i] < 0)
1414             {
1415               faces[i] = lowCellId;
1416               return;
1417             }
1418           if (faces[i] == lowCellId)
1419             return;
1420         }
1421     }
1422   ASSERT(0);
1423 }
1424
1425 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's
1426  * The pyramid is defined by the five points (0-4) where (0,1,2,3) is the base of the pyramid which,
1427  * using the right hand rule, forms a quadrilateral whose normal points in the direction of the
1428  * pyramid apex at vertex #4.
1429  * @see vtkPyramid.h in Filtering.
1430  * @param cellId volumeId in vtkUnstructuredGrid
1431  * @param facesWithNodes vector of face descriptors to be filled
1432  */
1433 void SMDS_DownPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1434 {
1435   // --- find point id's of the volume
1436
1437   vtkIdType npts = 0;
1438   vtkIdType *nodes; // will refer to the point id's of the volume
1439   _grid->GetCellPoints(cellId, npts, nodes);
1440
1441   // --- create all the ordered list of node id's for each face
1442
1443   facesWithNodes.nbElems = 5;
1444
1445   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1446   facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1447   facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1448   facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1449   facesWithNodes.elems[0].nbNodes = 4;
1450   facesWithNodes.elems[0].vtkType = VTK_QUAD;
1451
1452   facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1453   facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1454   facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1455   facesWithNodes.elems[1].nbNodes = 3;
1456   facesWithNodes.elems[1].vtkType = VTK_TRIANGLE;
1457
1458   facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1459   facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1460   facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1461   facesWithNodes.elems[2].nbNodes = 3;
1462   facesWithNodes.elems[2].vtkType = VTK_TRIANGLE;
1463
1464   facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1465   facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1466   facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1467   facesWithNodes.elems[3].nbNodes = 3;
1468   facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1469
1470   facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1471   facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1472   facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1473   facesWithNodes.elems[4].nbNodes = 3;
1474   facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1475 }
1476
1477 // ---------------------------------------------------------------------------
1478
1479 SMDS_DownQuadPyramid::SMDS_DownQuadPyramid(SMDS_UnstructuredGrid *grid) :
1480   SMDS_Down3D(grid, 5)
1481 {
1482   _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1483   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1484   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1485   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1486   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1487 }
1488
1489 SMDS_DownQuadPyramid::~SMDS_DownQuadPyramid()
1490 {
1491 }
1492
1493 void SMDS_DownQuadPyramid::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1494 {
1495 //   MESSAGE("SMDS_DownQuadPyramid::getOrderedNodesOfFace cellId = " << cellId);
1496   set<int> setNodes;
1497   setNodes.clear();
1498   for (int i = 0; i < orderedNodes.size(); i++)
1499     setNodes.insert(orderedNodes[i]);
1500   //MESSAGE("cellId = " << cellId);
1501
1502   vtkIdType npts = 0;
1503   vtkIdType *nodes; // will refer to the point id's of the volume
1504   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1505
1506   set<int> tofind;
1507   int ids[32] = { 0, 1, 2, 3, 5, 6, 7, 8,
1508                   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
1510   // Quadrangular face
1511   tofind.clear();
1512   for (int i = 0; i < 8; i++)
1513     tofind.insert(nodes[ids[i]]);
1514   if (setNodes == tofind)
1515     {
1516       for (int i = 0; i < 8; i++)
1517         orderedNodes[i] = nodes[ids[i]];
1518       return;
1519     }
1520   // Triangular faces
1521   for (int k = 0; k < 4; k++)
1522     {
1523       tofind.clear();
1524       for (int i = 0; i < 6; i++)
1525         tofind.insert(nodes[ids[8 + 6 * k + i]]);
1526       if (setNodes == tofind)
1527         {
1528           for (int i = 0; i < 6; i++)
1529             orderedNodes[i] = nodes[ids[8 + 6 * k + i]];
1530           return;
1531         }
1532     }
1533   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1534   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
1535   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1536 }
1537
1538 void SMDS_DownQuadPyramid::addDownCell(int cellId, int lowCellId, unsigned char aType)
1539 {
1540   //ASSERT((cellId >=0) && (cellId < _maxId));
1541   int *faces = &_cellIds[_nbDownCells * cellId];
1542   if (aType == VTK_QUADRATIC_QUAD)
1543     {
1544       if (faces[0] < 0)
1545         {
1546           faces[0] = lowCellId;
1547           return;
1548         }
1549       if (faces[0] == lowCellId)
1550         return;
1551     }
1552   else
1553     {
1554       //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1555       for (int i = 1; i < _nbDownCells; i++)
1556         {
1557           if (faces[i] < 0)
1558             {
1559               faces[i] = lowCellId;
1560               return;
1561             }
1562           if (faces[i] == lowCellId)
1563             return;
1564         }
1565     }
1566   ASSERT(0);
1567 }
1568
1569 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's
1570  * The ordering of the thirteen points defining the quadratic pyramid cell is point id's (0-4,5-12)
1571  * where point id's 0-4 are the five corner vertices of the pyramid; followed
1572  * by eight mid-edge nodes (5-12). Note that these mid-edge nodes lie on the edges defined by
1573  * 5(0,1), 6(1,2), 7(2,3), 8(3,0), 9(0,4), 10(1,4), 11(2,4), 12(3,4).
1574  * @see vtkQuadraticPyramid.h in Filtering.
1575  * @param cellId volumeId in vtkUnstructuredGrid
1576  * @param facesWithNodes vector of face descriptors to be filled
1577  */
1578 void SMDS_DownQuadPyramid::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1579 {
1580   // --- find point id's of the volume
1581
1582   vtkIdType npts = 0;
1583   vtkIdType *nodes; // will refer to the point id's of the volume
1584   _grid->GetCellPoints(cellId, npts, nodes);
1585
1586   // --- create all the ordered list of node id's for each face
1587
1588   facesWithNodes.nbElems = 5;
1589
1590   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1591   facesWithNodes.elems[0].nodeIds[1] = nodes[1];
1592   facesWithNodes.elems[0].nodeIds[2] = nodes[2];
1593   facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1594   facesWithNodes.elems[0].nodeIds[4] = nodes[5];
1595   facesWithNodes.elems[0].nodeIds[5] = nodes[6];
1596   facesWithNodes.elems[0].nodeIds[6] = nodes[7];
1597   facesWithNodes.elems[0].nodeIds[7] = nodes[8];
1598   facesWithNodes.elems[0].nbNodes = 8;
1599   facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1600
1601   facesWithNodes.elems[1].nodeIds[0] = nodes[0];
1602   facesWithNodes.elems[1].nodeIds[1] = nodes[1];
1603   facesWithNodes.elems[1].nodeIds[2] = nodes[4];
1604   facesWithNodes.elems[1].nodeIds[3] = nodes[5];
1605   facesWithNodes.elems[1].nodeIds[4] = nodes[10];
1606   facesWithNodes.elems[1].nodeIds[5] = nodes[9];
1607   facesWithNodes.elems[1].nbNodes = 6;
1608   facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
1609
1610   facesWithNodes.elems[2].nodeIds[0] = nodes[1];
1611   facesWithNodes.elems[2].nodeIds[1] = nodes[2];
1612   facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1613   facesWithNodes.elems[2].nodeIds[3] = nodes[6];
1614   facesWithNodes.elems[2].nodeIds[4] = nodes[11];
1615   facesWithNodes.elems[2].nodeIds[5] = nodes[10];
1616   facesWithNodes.elems[2].nbNodes = 6;
1617   facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
1618
1619   facesWithNodes.elems[3].nodeIds[0] = nodes[2];
1620   facesWithNodes.elems[3].nodeIds[1] = nodes[3];
1621   facesWithNodes.elems[3].nodeIds[2] = nodes[4];
1622   facesWithNodes.elems[3].nodeIds[3] = nodes[7];
1623   facesWithNodes.elems[3].nodeIds[4] = nodes[12];
1624   facesWithNodes.elems[3].nodeIds[5] = nodes[11];
1625   facesWithNodes.elems[3].nbNodes = 6;
1626   facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1627
1628   facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1629   facesWithNodes.elems[4].nodeIds[1] = nodes[0];
1630   facesWithNodes.elems[4].nodeIds[2] = nodes[4];
1631   facesWithNodes.elems[4].nodeIds[3] = nodes[8];
1632   facesWithNodes.elems[4].nodeIds[4] = nodes[9];
1633   facesWithNodes.elems[4].nodeIds[5] = nodes[12];
1634   facesWithNodes.elems[4].nbNodes = 6;
1635   facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1636 }
1637
1638 // ---------------------------------------------------------------------------
1639
1640 SMDS_DownPenta::SMDS_DownPenta(SMDS_UnstructuredGrid *grid) :
1641   SMDS_Down3D(grid, 5)
1642 {
1643   _cellTypes.push_back(VTK_QUAD);
1644   _cellTypes.push_back(VTK_QUAD);
1645   _cellTypes.push_back(VTK_QUAD);
1646   _cellTypes.push_back(VTK_TRIANGLE);
1647   _cellTypes.push_back(VTK_TRIANGLE);
1648 }
1649
1650 SMDS_DownPenta::~SMDS_DownPenta()
1651 {
1652 }
1653
1654 void SMDS_DownPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1655 {
1656   set<int> setNodes;
1657   setNodes.clear();
1658   for (int i = 0; i < orderedNodes.size(); i++)
1659     setNodes.insert(orderedNodes[i]);
1660   //MESSAGE("cellId = " << cellId);
1661
1662   vtkIdType npts = 0;
1663   vtkIdType *nodes; // will refer to the point id's of the volume
1664   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1665
1666   set<int> tofind;
1667 //int ids[18] = { 0, 2, 1,  3, 4, 5,   0, 1, 4, 3,   1, 2, 5, 4,   2, 0, 3, 5 };
1668   int ids[18] = { 0, 1, 2,  3, 5, 4,   0, 3, 4, 1,   1, 4, 5, 2,   2, 5, 3, 0 };
1669
1670   // Triangular faces
1671   for (int k = 0; k < 2; k++)
1672     {
1673       tofind.clear();
1674       for (int i = 0; i < 3; i++)
1675         tofind.insert(nodes[ids[3 * k + i]]);
1676       if (setNodes == tofind)
1677         {
1678           for (int i = 0; i < 3; i++)
1679             orderedNodes[i] = nodes[ids[3 * k + i]];
1680           return;
1681         }
1682     }
1683   // Quadrangular faces
1684   for (int k = 0; k < 3; k++)
1685     {
1686       tofind.clear();
1687       for (int i = 0; i < 4; i++)
1688         tofind.insert(nodes[ids[6 + 4 * k + i]]);
1689       if (setNodes == tofind)
1690         {
1691           for (int i = 0; i < 4; i++)
1692             orderedNodes[i] = nodes[ids[6 + 4 * k + i]];
1693           return;
1694         }
1695     }
1696   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1697   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1698   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1699 }
1700
1701 void SMDS_DownPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1702 {
1703   //ASSERT((cellId >=0) && (cellId < _maxId));
1704   int *faces = &_cellIds[_nbDownCells * cellId];
1705   if (aType == VTK_QUAD)
1706     for (int i = 0; i < 3; i++)
1707       {
1708         if (faces[i] < 0)
1709           {
1710             faces[i] = lowCellId;
1711             return;
1712           }
1713         if (faces[i] == lowCellId)
1714           return;
1715       }
1716   else
1717     {
1718       //ASSERT(aType == VTK_TRIANGLE);
1719       for (int i = 3; i < _nbDownCells; i++)
1720         {
1721           if (faces[i] < 0)
1722             {
1723               faces[i] = lowCellId;
1724               return;
1725             }
1726           if (faces[i] == lowCellId)
1727             return;
1728         }
1729     }
1730   ASSERT(0);
1731 }
1732
1733 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's.
1734  * A wedge or pentahedron consists of two triangular and three quadrilateral faces
1735  * and is defined by the six points (0-5) where (0,1,2) is the base of the wedge which,
1736  * using the right hand rule, forms a triangle whose normal points outward
1737  * (away from the triangular face (3,4,5)).
1738  * @see vtkWedge.h in Filtering
1739  * @param cellId volumeId in vtkUnstructuredGrid
1740  * @param facesWithNodes vector of face descriptors to be filled
1741  */
1742 void SMDS_DownPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1743 {
1744   // --- find point id's of the volume
1745
1746   vtkIdType npts = 0;
1747   vtkIdType *nodes; // will refer to the point id's of the volume
1748   _grid->GetCellPoints(cellId, npts, nodes);
1749
1750   // --- create all the ordered list of node id's for each face
1751
1752   facesWithNodes.nbElems = 5;
1753
1754   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1755   facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1756   facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1757   facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1758   facesWithNodes.elems[0].nbNodes = 4;
1759   facesWithNodes.elems[0].vtkType = VTK_QUAD;
1760
1761   facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1762   facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1763   facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1764   facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1765   facesWithNodes.elems[1].nbNodes = 4;
1766   facesWithNodes.elems[1].vtkType = VTK_QUAD;
1767
1768   facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1769   facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1770   facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1771   facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1772   facesWithNodes.elems[2].nbNodes = 4;
1773   facesWithNodes.elems[2].vtkType = VTK_QUAD;
1774
1775   facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1776   facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1777   facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1778   facesWithNodes.elems[3].nbNodes = 3;
1779   facesWithNodes.elems[3].vtkType = VTK_TRIANGLE;
1780
1781   facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1782   facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1783   facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1784   facesWithNodes.elems[4].nbNodes = 3;
1785   facesWithNodes.elems[4].vtkType = VTK_TRIANGLE;
1786 }
1787
1788 // ---------------------------------------------------------------------------
1789
1790 SMDS_DownQuadPenta::SMDS_DownQuadPenta(SMDS_UnstructuredGrid *grid) :
1791   SMDS_Down3D(grid, 5)
1792 {
1793   _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1794   _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1795   _cellTypes.push_back(VTK_QUADRATIC_QUAD);
1796   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1797   _cellTypes.push_back(VTK_QUADRATIC_TRIANGLE);
1798 }
1799
1800 SMDS_DownQuadPenta::~SMDS_DownQuadPenta()
1801 {
1802 }
1803
1804 void SMDS_DownQuadPenta::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1805 {
1806   set<int> setNodes;
1807   setNodes.clear();
1808   for (int i = 0; i < orderedNodes.size(); i++)
1809     setNodes.insert(orderedNodes[i]);
1810   //MESSAGE("cellId = " << cellId);
1811
1812   vtkIdType npts = 0;
1813   vtkIdType *nodes; // will refer to the point id's of the volume
1814   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1815
1816   set<int> tofind;
1817 //int ids[18] = { 0, 2, 1,  3, 4, 5,   0, 1, 4, 3,   1, 2, 5, 4,   2, 0, 3, 5 };
1818   int ids[36] = { 0, 1, 2, 6, 7, 8,  3, 5, 4, 11, 10, 9,
1819                   0, 3, 4, 1, 12, 9, 13, 6,   1, 4, 5, 2, 13, 10, 14, 7,   2, 5, 3, 0, 14, 11, 12, 8 };
1820
1821   // Triangular faces
1822   for (int k = 0; k < 2; k++)
1823     {
1824       tofind.clear();
1825       for (int i = 0; i < 6; i++)
1826         tofind.insert(nodes[ids[6 * k + i]]);
1827       if (setNodes == tofind)
1828         {
1829           for (int i = 0; i < 6; i++)
1830             orderedNodes[i] = nodes[ids[6 * k + i]];
1831           return;
1832         }
1833     }
1834   // Quadrangular faces
1835   for (int k = 0; k < 3; k++)
1836     {
1837       tofind.clear();
1838       for (int i = 0; i < 8; i++)
1839         tofind.insert(nodes[ids[12 + 8 * k + i]]);
1840       if (setNodes == tofind)
1841         {
1842           for (int i = 0; i < 8; i++)
1843             orderedNodes[i] = nodes[ids[12 + 8 * k + i]];
1844           return;
1845         }
1846     }
1847   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
1848   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
1849   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
1850 }
1851
1852 void SMDS_DownQuadPenta::addDownCell(int cellId, int lowCellId, unsigned char aType)
1853 {
1854   //ASSERT((cellId >=0) && (cellId < _maxId));
1855   int *faces = &_cellIds[_nbDownCells * cellId];
1856   if (aType == VTK_QUADRATIC_QUAD)
1857     for (int i = 0; i < 3; i++)
1858       {
1859         if (faces[i] < 0)
1860           {
1861             faces[i] = lowCellId;
1862             return;
1863           }
1864         if (faces[i] == lowCellId)
1865           return;
1866       }
1867   else
1868     {
1869       //ASSERT(aType == VTK_QUADRATIC_TRIANGLE);
1870       for (int i = 3; i < _nbDownCells; i++)
1871         {
1872           if (faces[i] < 0)
1873             {
1874               faces[i] = lowCellId;
1875               return;
1876             }
1877           if (faces[i] == lowCellId)
1878             return;
1879         }
1880     }
1881   ASSERT(0);
1882 }
1883
1884 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's
1885  * The quadratic wedge (or pentahedron) is defined by fifteen points.
1886  * The ordering of the fifteen points defining the cell is point id's (0-5,6-14)
1887  * where point id's 0-5 are the six corner vertices of the wedge, followed by
1888  * nine mid-edge nodes (6-14). Note that these mid-edge nodes lie on the edges defined by
1889  * (0,1), (1,2), (2,0), (3,4), (4,5), (5,3), (0,3), (1,4), (2,5).
1890  * @see vtkQuadraticWedge.h in Filtering
1891  * @param cellId volumeId in vtkUnstructuredGrid
1892  * @param facesWithNodes vector of face descriptors to be filled
1893  */
1894 void SMDS_DownQuadPenta::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
1895 {
1896   // --- find point id's of the volume
1897
1898   vtkIdType npts = 0;
1899   vtkIdType *nodes; // will refer to the point id's of the volume
1900   _grid->GetCellPoints(cellId, npts, nodes);
1901
1902   // --- create all the ordered list of node id's for each face
1903
1904   facesWithNodes.nbElems = 5;
1905
1906   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
1907   facesWithNodes.elems[0].nodeIds[1] = nodes[2];
1908   facesWithNodes.elems[0].nodeIds[2] = nodes[5];
1909   facesWithNodes.elems[0].nodeIds[3] = nodes[3];
1910   facesWithNodes.elems[0].nodeIds[4] = nodes[8];
1911   facesWithNodes.elems[0].nodeIds[5] = nodes[14];
1912   facesWithNodes.elems[0].nodeIds[6] = nodes[11];
1913   facesWithNodes.elems[0].nodeIds[7] = nodes[12];
1914   facesWithNodes.elems[0].nbNodes = 8;
1915   facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
1916
1917   facesWithNodes.elems[1].nodeIds[0] = nodes[1];
1918   facesWithNodes.elems[1].nodeIds[1] = nodes[2];
1919   facesWithNodes.elems[1].nodeIds[2] = nodes[5];
1920   facesWithNodes.elems[1].nodeIds[3] = nodes[4];
1921   facesWithNodes.elems[1].nodeIds[4] = nodes[7];
1922   facesWithNodes.elems[1].nodeIds[5] = nodes[14];
1923   facesWithNodes.elems[1].nodeIds[6] = nodes[10];
1924   facesWithNodes.elems[1].nodeIds[7] = nodes[13];
1925   facesWithNodes.elems[1].nbNodes = 8;
1926   facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
1927
1928   facesWithNodes.elems[2].nodeIds[0] = nodes[0];
1929   facesWithNodes.elems[2].nodeIds[1] = nodes[1];
1930   facesWithNodes.elems[2].nodeIds[2] = nodes[4];
1931   facesWithNodes.elems[2].nodeIds[3] = nodes[3];
1932   facesWithNodes.elems[2].nodeIds[4] = nodes[6];
1933   facesWithNodes.elems[2].nodeIds[5] = nodes[13];
1934   facesWithNodes.elems[2].nodeIds[6] = nodes[9];
1935   facesWithNodes.elems[2].nodeIds[7] = nodes[12];
1936   facesWithNodes.elems[2].nbNodes = 8;
1937   facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
1938
1939   facesWithNodes.elems[3].nodeIds[0] = nodes[0];
1940   facesWithNodes.elems[3].nodeIds[1] = nodes[1];
1941   facesWithNodes.elems[3].nodeIds[2] = nodes[2];
1942   facesWithNodes.elems[3].nodeIds[3] = nodes[6];
1943   facesWithNodes.elems[3].nodeIds[4] = nodes[7];
1944   facesWithNodes.elems[3].nodeIds[5] = nodes[8];
1945   facesWithNodes.elems[3].nbNodes = 6;
1946   facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
1947
1948   facesWithNodes.elems[4].nodeIds[0] = nodes[3];
1949   facesWithNodes.elems[4].nodeIds[1] = nodes[4];
1950   facesWithNodes.elems[4].nodeIds[2] = nodes[5];
1951   facesWithNodes.elems[4].nodeIds[3] = nodes[9];
1952   facesWithNodes.elems[4].nodeIds[4] = nodes[10];
1953   facesWithNodes.elems[4].nodeIds[5] = nodes[11];
1954   facesWithNodes.elems[4].nbNodes = 6;
1955   facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_TRIANGLE;
1956 }
1957
1958 // ---------------------------------------------------------------------------
1959
1960 SMDS_DownHexa::SMDS_DownHexa(SMDS_UnstructuredGrid *grid) :
1961   SMDS_Down3D(grid, 6)
1962 {
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);
1968   _cellTypes.push_back(VTK_QUAD);
1969 }
1970
1971 SMDS_DownHexa::~SMDS_DownHexa()
1972 {
1973 }
1974
1975 void SMDS_DownHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
1976 {
1977   set<int> setNodes;
1978   setNodes.clear();
1979   for (int i = 0; i < orderedNodes.size(); i++)
1980     setNodes.insert(orderedNodes[i]);
1981   //MESSAGE("cellId = " << cellId);
1982
1983   vtkIdType npts = 0;
1984   vtkIdType *nodes; // will refer to the point id's of the volume
1985   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
1986
1987   set<int> tofind;
1988 //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};
1989   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};
1990   for (int k = 0; k < 6; k++) // loop on the 6 faces
1991     {
1992       tofind.clear();
1993       for (int i = 0; i < 4; i++)
1994         tofind.insert(nodes[ids[4 * k + i]]); // node ids of the face i
1995       if (setNodes == tofind)
1996         {
1997           for (int i = 0; i < 4; i++)
1998             orderedNodes[i] = nodes[ids[4 * k + i]];
1999           return;
2000         }
2001     }
2002   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
2003   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2004   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2005   MESSAGE(nodes[4] << " " << nodes[5] << " " << nodes[6] << " " << nodes[7]);
2006 }
2007
2008 void SMDS_DownHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
2009 {
2010   //ASSERT((cellId >=0)&& (cellId < _maxId));
2011   int *faces = &_cellIds[_nbDownCells * cellId];
2012   for (int i = 0; i < _nbDownCells; i++)
2013     {
2014       if (faces[i] < 0)
2015         {
2016           faces[i] = lowCellId;
2017           return;
2018         }
2019       if (faces[i] == lowCellId)
2020         return;
2021     }
2022   ASSERT(0);
2023   // MESSAGE("-------------------------------------> trop de faces ! " << cellId << " " << lowCellId);
2024 }
2025
2026 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's
2027  * The hexahedron is defined by the eight points (0-7), where (0,1,2,3) is the base
2028  * of the hexahedron which, using the right hand rule, forms a quadrilateral whose normal
2029  * points in the direction of the opposite face (4,5,6,7).
2030  * @see vtkHexahedron.h in Filtering
2031  * @param cellId volumeId in vtkUnstructuredGrid
2032  * @param facesWithNodes vector of face descriptors to be filled
2033  */
2034 void SMDS_DownHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2035 {
2036   // --- find point id's of the volume
2037
2038   vtkIdType npts = 0;
2039   vtkIdType *nodes; // will refer to the point id's of the volume
2040   _grid->GetCellPoints(cellId, npts, nodes);
2041
2042   // --- create all the ordered list of node id's for each face
2043
2044   facesWithNodes.nbElems = 6;
2045
2046   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2047   facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2048   facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2049   facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2050   facesWithNodes.elems[0].nbNodes = 4;
2051   facesWithNodes.elems[0].vtkType = VTK_QUAD;
2052
2053   facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2054   facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2055   facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2056   facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2057   facesWithNodes.elems[1].nbNodes = 4;
2058   facesWithNodes.elems[1].vtkType = VTK_QUAD;
2059
2060   facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2061   facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2062   facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2063   facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2064   facesWithNodes.elems[2].nbNodes = 4;
2065   facesWithNodes.elems[2].vtkType = VTK_QUAD;
2066
2067   facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2068   facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2069   facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2070   facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2071   facesWithNodes.elems[3].nbNodes = 4;
2072   facesWithNodes.elems[3].vtkType = VTK_QUAD;
2073
2074   facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2075   facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2076   facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2077   facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2078   facesWithNodes.elems[4].nbNodes = 4;
2079   facesWithNodes.elems[4].vtkType = VTK_QUAD;
2080
2081   facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2082   facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2083   facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2084   facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2085   facesWithNodes.elems[5].nbNodes = 4;
2086   facesWithNodes.elems[5].vtkType = VTK_QUAD;
2087 }
2088
2089 // ---------------------------------------------------------------------------
2090
2091 SMDS_DownQuadHexa::SMDS_DownQuadHexa(SMDS_UnstructuredGrid *grid) :
2092   SMDS_Down3D(grid, 6)
2093 {
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);
2099   _cellTypes.push_back(VTK_QUADRATIC_QUAD);
2100 }
2101
2102 SMDS_DownQuadHexa::~SMDS_DownQuadHexa()
2103 {
2104 }
2105
2106 void SMDS_DownQuadHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
2107 {
2108   set<int> setNodes;
2109   setNodes.clear();
2110   for (int i = 0; i < orderedNodes.size(); i++)
2111     setNodes.insert(orderedNodes[i]);
2112   //MESSAGE("cellId = " << cellId);
2113
2114   vtkIdType npts = 0;
2115   vtkIdType *nodes; // will refer to the point id's of the volume
2116   _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
2117
2118   set<int> tofind;
2119   //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};
2120   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,
2121                   4, 0, 1, 5,16, 8,17,12,   5, 1, 2, 6,17, 9,18,13,   6, 2, 3, 7,18,10,19,14};
2122   for (int k = 0; k < 6; k++)
2123     {
2124       tofind.clear();
2125       for (int i = 0; i < 8; i++)
2126         tofind.insert(nodes[ids[8 * k + i]]);
2127       if (setNodes == tofind)
2128         {
2129           for (int i = 0; i < 8; i++)
2130             orderedNodes[i] = nodes[ids[8 * k + i]];
2131           return;
2132         }
2133     }
2134   MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
2135   MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
2136   MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
2137 }
2138
2139 void SMDS_DownQuadHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)
2140 {
2141   //ASSERT((cellId >=0)&& (cellId < _maxId));
2142   int *faces = &_cellIds[_nbDownCells * cellId];
2143   for (int i = 0; i < _nbDownCells; i++)
2144     {
2145       if (faces[i] < 0)
2146         {
2147           faces[i] = lowCellId;
2148           return;
2149         }
2150       if (faces[i] == lowCellId)
2151         return;
2152     }
2153   ASSERT(0);
2154 }
2155
2156 /*! Create a list of faces described by a vtk Type and  an ordered set of Node Id's
2157  * The ordering of the twenty points defining the quadratic hexahedron cell is point id's (0-7,8-19)
2158  * where point id's 0-7 are the eight corner vertices of the cube, followed by twelve mid-edge nodes (8-19).
2159  * Note that these mid-edge nodes lie on the edges defined by
2160  * (0,1), (1,2), (2,3), (3,0), (4,5), (5,6), (6,7), (7,4), (0,4), (1,5), (2,6), (3,7).
2161  * @see vtkQuadraticHexahedron.h in Filtering
2162  * @param cellId volumeId in vtkUnstructuredGrid
2163  * @param facesWithNodes vector of face descriptors to be filled
2164  */
2165 void SMDS_DownQuadHexa::computeFacesWithNodes(int cellId, ListElemByNodesType& facesWithNodes)
2166 {
2167   // --- find point id's of the volume
2168
2169   vtkIdType npts = 0;
2170   vtkIdType *nodes; // will refer to the point id's of the volume
2171   _grid->GetCellPoints(cellId, npts, nodes);
2172
2173   // --- create all the ordered list of node id's for each face
2174
2175   facesWithNodes.nbElems = 6;
2176
2177   facesWithNodes.elems[0].nodeIds[0] = nodes[0];
2178   facesWithNodes.elems[0].nodeIds[1] = nodes[1];
2179   facesWithNodes.elems[0].nodeIds[2] = nodes[2];
2180   facesWithNodes.elems[0].nodeIds[3] = nodes[3];
2181   facesWithNodes.elems[0].nodeIds[4] = nodes[8];
2182   facesWithNodes.elems[0].nodeIds[5] = nodes[9];
2183   facesWithNodes.elems[0].nodeIds[6] = nodes[10];
2184   facesWithNodes.elems[0].nodeIds[7] = nodes[11];
2185   facesWithNodes.elems[0].nbNodes = 8;
2186   facesWithNodes.elems[0].vtkType = VTK_QUADRATIC_QUAD;
2187
2188   facesWithNodes.elems[1].nodeIds[0] = nodes[4];
2189   facesWithNodes.elems[1].nodeIds[1] = nodes[5];
2190   facesWithNodes.elems[1].nodeIds[2] = nodes[6];
2191   facesWithNodes.elems[1].nodeIds[3] = nodes[7];
2192   facesWithNodes.elems[1].nodeIds[4] = nodes[12];
2193   facesWithNodes.elems[1].nodeIds[5] = nodes[13];
2194   facesWithNodes.elems[1].nodeIds[6] = nodes[14];
2195   facesWithNodes.elems[1].nodeIds[7] = nodes[15];
2196   facesWithNodes.elems[1].nbNodes = 8;
2197   facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_QUAD;
2198
2199   facesWithNodes.elems[2].nodeIds[0] = nodes[0];
2200   facesWithNodes.elems[2].nodeIds[1] = nodes[1];
2201   facesWithNodes.elems[2].nodeIds[2] = nodes[5];
2202   facesWithNodes.elems[2].nodeIds[3] = nodes[4];
2203   facesWithNodes.elems[2].nodeIds[4] = nodes[8];
2204   facesWithNodes.elems[2].nodeIds[5] = nodes[17];
2205   facesWithNodes.elems[2].nodeIds[6] = nodes[12];
2206   facesWithNodes.elems[2].nodeIds[7] = nodes[16];
2207   facesWithNodes.elems[2].nbNodes = 8;
2208   facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_QUAD;
2209
2210   facesWithNodes.elems[3].nodeIds[0] = nodes[1];
2211   facesWithNodes.elems[3].nodeIds[1] = nodes[2];
2212   facesWithNodes.elems[3].nodeIds[2] = nodes[6];
2213   facesWithNodes.elems[3].nodeIds[3] = nodes[5];
2214   facesWithNodes.elems[3].nodeIds[4] = nodes[9];
2215   facesWithNodes.elems[3].nodeIds[5] = nodes[18];
2216   facesWithNodes.elems[3].nodeIds[6] = nodes[13];
2217   facesWithNodes.elems[3].nodeIds[7] = nodes[17];
2218   facesWithNodes.elems[3].nbNodes = 8;
2219   facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_QUAD;
2220
2221   facesWithNodes.elems[4].nodeIds[0] = nodes[2];
2222   facesWithNodes.elems[4].nodeIds[1] = nodes[6];
2223   facesWithNodes.elems[4].nodeIds[2] = nodes[7];
2224   facesWithNodes.elems[4].nodeIds[3] = nodes[3];
2225   facesWithNodes.elems[4].nodeIds[4] = nodes[18];
2226   facesWithNodes.elems[4].nodeIds[5] = nodes[14];
2227   facesWithNodes.elems[4].nodeIds[6] = nodes[19];
2228   facesWithNodes.elems[4].nodeIds[7] = nodes[10];
2229   facesWithNodes.elems[4].nbNodes = 8;
2230   facesWithNodes.elems[4].vtkType = VTK_QUADRATIC_QUAD;
2231
2232   facesWithNodes.elems[5].nodeIds[0] = nodes[3];
2233   facesWithNodes.elems[5].nodeIds[1] = nodes[7];
2234   facesWithNodes.elems[5].nodeIds[2] = nodes[4];
2235   facesWithNodes.elems[5].nodeIds[3] = nodes[0];
2236   facesWithNodes.elems[5].nodeIds[4] = nodes[19];
2237   facesWithNodes.elems[5].nodeIds[5] = nodes[15];
2238   facesWithNodes.elems[5].nodeIds[6] = nodes[16];
2239   facesWithNodes.elems[5].nodeIds[7] = nodes[11];
2240   facesWithNodes.elems[5].nbNodes = 8;
2241   facesWithNodes.elems[5].vtkType = VTK_QUADRATIC_QUAD;
2242 }
2243
2244 // ---------------------------------------------------------------------------
2245